PEIiMIT  FlUY  lEGIIIlf  PRIHHICIION 


UNITED  STATES  AIR  FORCE 
AIR  UNIVERSITY 


AIR  FORCE  INSTITUTE  OF  TECHNOLOGY 

Wriglit-Fott«rton  Air  Fere*  Ras«,Oklo 


GNE/PH/76D-1 


THESIS 


aUE/PH/76D-l 


Gene  P.  Bender 
Lt  USN 


Approved  for  public  release;  distribution  unlimited 


I 


INTERACTION  OF  RELATIVISTIC 
PROTONS  WITH  MATTER 


gNE/PH/76D-l  I 


INTERACTION  OF  RELATIVISTIC 
PROTONS  WITH  MATTER  . 


Presented  to  the  Faculty  of  the  School  of  Engineering 


of  the  Air  Force  Institute  of  Technology 


Air  University 


in  Partial  Fulfillment  of  the 


Requirements  for  the  Degree  of 


Master  of  Science 


hr 


UTSS 

CK 

atti'erut 


ef 

SKi^w  Q 


*TlfiC4TXil  i 




Gene  P. /Bender 


\ B.S. 


GraRu^te  Engineering  Physics 
n Deceator  39^ 


Approved  for  public  release;  distribution  unlimited. 


lited. 


Preface 


In  this  thesis,  an  attempt  was  made  to  use  a simple 
and  direct  approach  to  solve  a complicated  physical  problem. 
As  with  most  problems  of  this  type,  the  analysis  raised  more 
questions  than  it  answered.  I feel  that  the  problem  of  col- 
lisionless  Interactions  deserves  further  detailed  study,  in 
order  to  truly  establish  the  mechanism  as  an  important  con- 
tributor to  energy  deposition  in  a target. 

I am  deeply  grateful  to  Capt.  P.  Nielsen  and  to  Dr. 

D.  Shankland  for  both  their  guidance  and  instruction  dviring 
the  course  of  this  work  and  to  those  at  home  who  are  patient 
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Abstract 

The  mechanisms  associated  with  the  interaction  of  a 
relativistic  proton  beam  with  a target  material  are  sum- 
marized, and  a simple  analytic  calculation  for  each  mech- 
anism relates  its  importance  in  the  deposition  of  energy. 
Because  it  was  fovind  that  possible  significant  contributions 
coiild  be  related  to  collective  effects,  the  first  three 
moments  of  the  Boltzmann-Vlasov  equation  and  the  equations 
necessary  to  describe  the  internal  and  external  electromag- 
netic fields  associated  with  the  target  are  developed. 

These  equations  are  then  formulated  into  a computer  program 
to  describe  the  collective  interaction  process.  Sample 
r^lns  demonstrate  the  ability  of  the  computer  model  to  pre- 
dict plasma  oscillations  at  the  proper  frequency.  Al- 
though an  instability  developed  in  the  runs  with  the  beam 
approaching  the  target,  the  preliminary  analysis  was  ac- 
complished, and  the  need  for  future  detailed  analysis  of 
this  process  was  established. 


I . Introduction 


In  the  past  several  years,  there  has  been  considerable 
interest  in  the  interaction  of  relativistic  electrons  and 
ions  with  matter.  This  interest  has  been  stirred  in  part 

i 

by  the  prospect  of  major  breakthroughs  in  the  ability  to  | 

I 

accelerate  charged  particles  to  GeV  energies  with  accelera-  | 

tors  of  reasonable  size.  Although  most  of  the  research  be- 
ing done  today  concerns  construction  of  the  accelerators  or 
the  propagation  of  the  beam  in  the  accelerator  and  through 
the  atmosphere,  some  work  is  also  being  done  on  the  mechan- 
ism of  interaction  of  these  beams  with  the  target. 

The  known  interaction  mechanisms  associated  with  the 
charged  particle  beam  striking  the  target  can  be  divided 
into  the  following  two  major  categories: 

1.  Collislonal  interactions 

2.  Collective  (sometimes  termed  collisionless) 

Interactions 

Collislonal  interactions  include  ionization  losses  due 
to  inelalastic  collisions  with  the  atomic  electrons,  elastic 
interactions  with  the  nucleus,  Bremsstrahlung  produced  as 
the  result  of  collisions,  and  inelastic  nuclear  interactions. 

Collective  interactions  are  those  associated  with  the 
system  as  a group  and  not  individual  elements.  External 
forces  can  induce  oscillations  within  a plasma,  and  these 
oscillations  may  then  be  converted  into  random  motion  or 
heat  through  instabilities  such  as  two-stream  or  return-cur- 
rent heating. 

To  date,  considerable  research  has  been  done  on  the  j 
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collisional  interactions  of  protons  and  electrons  and  col- 
lective effects  of  electron  beams;  however,  very  little  has 
been  done  on  the  collective  interactions  of  protom  beams 
with  matter.  In  this  thesis,  I propose,  first,  to  svirvey 
the  major  interaction  mechanisms  associated  with  the  proton 
beam;  secondly,  to  develop  the  equations  necessary  to 
model  the  collective  interaction  mechanism;  and  thirdly, 
to  use  these  equations  to  investigate  the  physical  aspects 
of  the  relativistic  proton  beam  passing  through  a target 
material. 

History 

In  the  early  1950's,  the  charged  particle  beam  was  in- 
vestigated extensively  as  to  the  prospects  of  using  the  beam 
for  research  in  nuclear  structure  by  bombardment  of  nuclei 
with  protons,  investigation  of  high  energy  interaction 
mechanisms,  transmission  of  high  energy  charged  particles 
through  waveguides,  and  also  as  a weapon  system.  After 
considerable  time  and  effort,  the  state  of  the  art  limi- 
tations on  accelerators  and  electronic  compoents  slowed 
the  effort  to  a trickle.  Today,  with  a new  class  of  ac- 
celerators on  the  production  horizon,  the  investigation 
of  charged  ion  beams  of  relativistic  energies  has  begun 
with  renewed  interest. 

At  present,  numerous  institutions  are  undertaking  the 
development  of  collective  ion  accelerators  which  will  be 
capable  of  producing  the  energies  and  quantities  of  protons 
or  heavier  ions  necessary  for  propagation  through  the  atmos- 
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phere . 

Not  only  institutions  within  the  United  States,  but 
the  Soviet  Union  is  also  actively  engaged  in  studying  the 
interaction  and  propagation  mechanisms  of  the  charged  par- 
ticle beam.  The  topic  of  greatest  concern  in  much  of 
today's  literature  involves  the  choice  of  the  electron 
or  the  proton  for  most  desirable  propagational  and  in- 
teraction properties.  There  are  pros  and  cons  on  both 
sides  of  the  question;  however,  the  consensus  is,  if  the 
proton  beam  can  be  produced,  then  it  is  the  superior  par- 
ticle for  Interaction.  Thus,  at  this  point,  the  dis- 
cussion will  be  limited  to  the  proton  beam  except  where 
noted . 

Analysis 

Assumption.  The  most  basic  assumption  in  this  analysis 
is  the  beam  itself.  The  assumption  is  that  a proton  accel- 
erator capable  of  providing  high  enough  energies  and  den- 
sities to  propagate  a beam  over  reasonable  distance  can  be 
produced. 

The  second  major  assumption  concerns  the  ability  of 
n\americal  solutions  to  model  accurately  the  real  world  pro- 
cesses involved  in  the  interaction.  Some  assumptions  will 
be  stated  later  which  have  been  chosen  to  permit  a simpli- 
fication of  the  model  and  equations  used  in  the  ntamerical 
calc\ilations . These  assumptions  may  produce  slight  de- 
riations  from  the  real  world  problem,  but  the  overall 
trend  and  physics  should  still  be  preserved. 
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Facts . The  premise  used  for  the  modeling  of  the  col- 
lective and  collisional  interaction  mechanisms  has  come  from 
experimental  studies  done  on  very  small  scales.  Although 
^hese  small-scale  experiments  may  not  lead  to  a direct  con- 
version to  experiments  on  larger  scales,  the  physics  of  the 
interactions  should  remain  the  same. 

Electron  beam  propagation  and  Interaction  processes 
have  been  studied  extensively  over  the  past  twenty  years 
under  various  research  topics.  Even  though  single  proton 
beam  interactions  are  well  categornized,  limited  amovints  of 
research  have  been  done  with  relativistic  proton  beams 
with  large  number  densities  in  the  beam.  In  the  past  fif- 
teen years,  the  Soviet  Union  has  devoted  large  amounts  of 
time  and  effort  in  studying  the  relativistic  proton  beam  in- 
action mechanism,  and  is  today  the  probable  leader  in 
understanding  the  mechanisms  of  collisional  processes. 

Collective  Interactions  have  been  used  to  explain  such 
things  as  the  effects  on  communication  systems  by  solar 
flares  for  the  past  thirty  years.  The  basic  theory  of 
collective  interactions  is  well  documented;  however,  appli- 
cation of  this  phenomenon  as  a proton  beam  interaction 
mechanism  has  not  been  attempted. 

Sutiimarv 

This  thesis  will  summarize  the  various  types  of  inter- 
actions associated  with  a proton  beam  and  look  in  detail  at 
the  process  of  collective  interactions.  The  results  will 
tend  to  support  the  collective  Interaction  mechanism  as 


an  important  part  of  the  overall  interaction  process  or 
eliminate  the  mechanism  as  being  a minor  contributor  to  the 
deposition  of  energy  in  a target  material. 


A particle  passing  through  matter  makes  collisions  with 


both  the  atomic  electrons  and  the  nuclei  of  the  target  ma- 
terial. By  collision,  it  is  meant  that  the  very  strong 
short-range  forces  associated  with  interactions  come  into 
play.  Short-range  forces  involve  the  Interaction  of  either 
the  nuclear  forces  or  the  coulomblc  fields  associated  with 
the  charges  of  the  two  colliding  particles.  The  range  of 
these  collisions  is  limited  to  several  atomic  radii.  The 
proton,  being  some  two  thousand  times  heavier  than  the 
electron,  will  Impart  considerable  energy  to  the  electron, 
through  the  coulomb  field  collision,  knocking  it  from  orbit 
and  leaving  an  ion  behind.  This  process,  known  as  ioniza- 
tion, is  a major  contributor  to  degradation  of  the  beam 
energy.  Collisions  with  the  nucleus  can  take  place  via 
either  elastic  or  inelastic  mechanisms.  In  elastic  col- 
lisions of  a proton  with  a nucleus,  the  total  kinetic 
energy  of  the  system  remains  the  same  before  and  after  the 
collision;  thus,  no  energy  is  stored  in  the  nucleus,  and 
the  collision  is  termed  elastic.  If  the  number  of  par^ 
tides  Involved  and  their  energy  differ  after  collision 
or  the  total  kinetic  energy  is  not  conserved,  the  col- 
lision is  termed  inelastic.  Bremsstrahlung  is  produced  as 
the  result  of  acceleration  or  deceleration  of  charged 
particles  as  the  result  of  collisions.  This  radiation 
will  be  stopped  in  the  target  material  producing  heat. 

The  following  discussion  considers  each  of  these  mechan- 


6 


RELATIVE  IONIZATION 


isms  in  more  detail  and  tries  to  bring  them  into  per- 
spective as  to  importance  in  energy  deposition  in  the  target. 


:| 

Ionization  Losses 

ii  

I 

! V/hen  heavy  particles,  such  as  protons  with  high  energy, 

' pass  through  matter,  they  do  so  in  relatively  straight  lines 

f|  dissipating  a portion  of  their  energy  through  a series  of 

collision  with  the  atomic  electrons.  This  energy  dissipa- 
tion does  not  take  place  in  a linear  fashion,  but  as  can  be 
seen  by  Fig.  1,  depends  upon  the  velocity  of  the  particle. 

^ This  figure  relates  the  number  of  ion  pairs  per  millimeter 

of  path  length  (specific  ionization)  formed  as  a result 


of  collisions  of  the  incident  particle  with  the  orbital 
electrons  in  the  target  material. 


r 

r 
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It  can  be  seen  that  more  ion  pairs  are  formed  and  thus  the 
largest  energy  losses  occur  at  the  end  of  the  particle's 
range. 

Another  term  frequently  referred  to  in  the  literature 
when  talking  about  ionization  losses  is  "stopping  power", 

S i which  is  defined  as  the  energy  lost  by  a particle 
per  unit  length  of  travel  in  the  stopping  substance  (Ref. 

2) 

S (E)  = - (1) 

\ 

where  E is  the  kinetic  energy  of  the  particle,  and  could  be 
related  to  the  slope  of  Fig  1.  Stopping  power,  though 
usually  thought  of  as  a property  of  the  medium,  relates 
the  energy  transferred  in  a collision  of  a particle  of 
charge  (?e)  and  mass  (A^)  with  the  orbital  electrons  of 
the  stopping  substance.  Theoretically,  no  nuclear  colli- 
sions are  taken  into  account  in  this  term. 

If  the  velocity  of  the  Incident  particle  is  large  in 
comparison  with  the  orbital  velocity  of  the  electrons j 

5 

(electrons  at  room  temperature  have  a velocity  of  ^1  x 10 
m/sHhe  electrons  can  be  considered  at  rest  (Ref.  3)*  If 
we  assxime  that  the  momentum  transfer,  AP,  is  s\ifficiently 
small  that  the  incident  particle  is  not  deflected  from  its 
straight-line  path  and  that  the  recoiling  electron  does 
not  move  appreciably  during  the  collision,  then  id P is  a 
function  of  the  electric  field  of  the  incident  particle 
at  the  position  of  the  electron. 
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Magnetic  effects  are  negligible  if  the  electron  is  essen- 
tially at  rest. 


lAf 


V = rf/eioexit 

Puf-Ti'ttt. 

E r A|  c*  r 

M<  « 
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Only  the  transverse  component  of  the  electric  field  has  a 
non-vanishing  time  integral  in  one  dimension  and  is  there- 


fore given  by: 


Yg  b 


Where  and  x j ( i - v^/c*-) (Ref  3)* 

Thus  the  momentum  impulse,  also  only  in  the  transverse 
direction,  has  magnitude 

AP=  j (3) 

The  energy  transferred  to  the  electron  is 


7"  = 1a^ 


Mtb'-Y’- 


From  Eq.  (4)  we  see  that  the  energy  transferred, ^ D 
depends  upon  the  charge  on  the  incident  particle,  j.  in- 
versely as  the  square  of  the  impact  parameter,  b , and 
is  independent  of  the  incident  particle. 

A very  fast  particle  passing  through  matter  “sees” 
electrons  at  various  impact  parameters . The  number  of 
possible  collision  sites  ( diY)  , located  at  impact  para- 
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meters  between  (b)  and  (b  + db),  in  a thickness  dx  of  the 


material  is 


dN  - N a dK 


where  i is  the  atomic  nvimber  of  the  stopping  material. 

To  find  the  energy  lost  per  unit  distance  traveled  by 
the  incident  particle,  multiply  Eq.  (5)  by  the  energy  trans- 
ferred AT(b)  and  integrate  over  the  impact  parameters. 

Ni  2TT  ‘(t 


dx  ■ Ale  \/* 


or. 


fa 


b « 


The  correct  choice  of  b^,,^  and  b^.,*  must  be  based  on 
quantum  mechanical  considerations  as  was  done  by  Bethe 
(Ref.  20),  but  a rough  indication  of  the  result  can  be 
derived  classically.  The  lower  limit ,' i>/vnA/  , insures 
from  Eq.  (h)  and  on  noting  that  the  maximum  momentum  trans- 
ferred to  an  electron  in  a classical  central  collision  is 
*2  V , that  . 

The  upper  limit,  , is  determined  by  the  finite 
binding  energies  of  the  electrons  in  an  atom.  A bound 
electron  cannot  accept  arbitrarily  small  amounts  of  energy 
but  must  absorb  sufficient  energy  to  be  raised  to  an  un- 
filled energy  state.  We  may  introduce  effective  or  average 
excitation  energy,  I,  averaged  over  all  the  electrons  in 
the  atoms  of  the  absorber  and  contend  that  there  is  no 
energy  transferred  \inless  X , i.  e.; 

2 ^’•e"  j X (7) 

Al.  ■ ' 

toeref ore , 


J,  i b„.. (*) 

Substitution  of  these  limits  into  Eq  (6)  yields  the  approx- 
imate expression 


which  differs  from  the  corresponding  quantxim  mechanical  ex- 
pression only  in  the  occurrence  of  the  power  i in  the  log- 
arithmic term.  The  quantum  mechanical  expression  for  stop- 
ping power  is  (Ref.  2) 

- ^ A'-'  ^ 

Empirical  values  for  I,  the  average  excitation  potential,  are 
given  by 

Xsi  9.1Z  (1  -h  ) ev.  (/') 

Sternhelmer  (Ref  6)  corrected  the  Bethe  equation  [|Eq. 
(,10)J  to  Include  both  relativistic  particles  and  a better 
fit  to  experimental  data  at  lower  energies. 

-4f--  ^ <-> 

where  ^ = % • 

P = a.  correction  for  non-participation  of  the  K,  1, 

. . . shell  electrons  if  the  transferred  energy  of 
the  incident  particle  is  below  that  required  to 
overcome  the  nuclear  forces. 

6 = a correction  for  polarization  and  conductivity 
in  dense  materials.  This  factor  takes  into  ac- 
count the  reduction  of  ionization  losses  at  high 
energies  due  to  polarization  of  the  medixun  by  the 
effective  field  of  the  passing  charged  particles. 
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This  effect  becomes  important  above  kinetic 
energies  greater  than  2 GeV.  The  value  of 
§ for  specific  materials  and  energies  of  inci- 
dent particles  is  available  for  most  materials 
(Ref  5 and  6) . 

Elastic  Scattering 

In  the  interaction  process,  elastic  scattering  is  a 
minor  contributor  to  the  deposition  of  energy  in  the  tar- 
get material.  However,  the  transverse  spreading  of  the 
beam  during  propagation  of  the  protons  to  the  target  is 
of  significance,  thus  a brief  summaryof  the  process  will 
be  given. 

Elastic  scattering  requires  the  total  kinetic  energy 
of  the  system,  projectile  and  target,  to  be  the  same  before 
the  collision  as  after.  Some  kinetic  energy  is  usually 
transferred  from  the  Incident  particle  to  the  target  nu- 
cleus, but  the  latter  is  left  in  the  same  internal  or  nu- 
clear state,  as  before  the  collision.  A collision  between 
two  billiard  balls  , in  which  neither  is  damaged  or  changed 
otherwise,  is  an  elastic  collision. 

In  the  propagation  of  the  beam  through  a medium,  the 
elastic  collisions  will  impart  transverse  momentum  to  the 
scattered  proton.  This  transverse  momentum  will  tend  to 
spread  the  beam  in  radius  and  therefore  some  particles  will 
not  strike  the  target.  Once  the  beam  has  reached  the  target, 
the  elastic  collisions  within  the  target  will  involve  a 
spreading  of  the  deposition  energy.  This  spreading  will 
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tend  to  decrease  local  heating  at  the  beam  center  and  thus 
reduce  the  damage  caused  by  thermal  heating.  The  cross 
sections  for  elastic  interactions  are  slightly  less  than 
300  millibarns  (Ref.  7)  for  particles  in  the  GeV  energy 
range  in  aluminum.  Thus,  elastic  collision  will  not  re- 
sult in  significant  beam  spreading  within  the  target 
material.  The  equations  which  describe  this  spread  can  be 
fotind  in  Refs.  3 and  7>  but  under  the  assumption  that  the 
beam  has  already  reached  the  target,  this  phenomenon  will 
be  left  to  those  studying  propagational  processes. 

Nuclear  Interactions 

At  proton  energies  of  several  h\mdred  MeV  and  greater, 
the  collision  of  protons  with  the  nucleus  involves  a strong 
interaction  of  the  incident  particles  with  the  nucleons. 

This  interaction  leads  to  the  production  of  cascade  par- 
ticles, including  pi-mesons.  As  the  proton  energy  in- 
creases further  (near  6 GeV),  the  production  of  antipar- 
ticles also  becomes  energetically  possible.  The  cascade 
particles  result  from  the  plmary  collisions  near  the  edge  of 
the  nucleus  where  the  primary  particles  enter  , and  each 
secondary  particle  also  has  the  capability  of  causing 
additional  interactions  of  sufficient  energy  is  present. 

The  cross  sections  for  these  inelastic  collisions 
of  relatlvisltc  protons  are  approximately  constant  in  the 
10  GeV  range,  as  shown  by  Fig.  2. 
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Steenberg  (Ref  21)  has  developed  an  empirical  for- 
mula for  these  high  energy  proton  cross  sections  as  a 
fiinction  of  atomic  number  of  the  target  materials. 

where  = 38  millibarns,  and  A is  the  atomic 

number  of  the  target  material.  With  such  small  cross 
sections,  either  the  stopping  material  will  have  to  be 
greater  than  normal  density,  or  the  nvimber  of  protons  in 
the  beam  will  have  to  be  exceedingly  large  in  order  to 
• deposit  an  appreciable  amount  of  energy  within  a reason- 

able range. 

The  interaction  of  high  energy  protons  with  the  nucle- 
j ons  may  produce  a variety  of  particle  species  and  re- 

^ action  mechanisms  as  shown  by  Fig.  3*  This  figure  pictor- 

ially  describes  the  cascading  process  in  which  second- 

14 
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ary  particles  interact  with  surrounding  nucleons  and  pro- 
duce mviltiple  interactions.  The  study  of  the  nuclear  cas- 
cade process  is  the  subject  of  considerable  amounts  of  pre- 
sent-day research  in  both  the  United  States  and  the  Soviet 
Union.  The  problems  associated  with  the  modeling  of  these 
cascades  have  produced  extensive  discussion  as  to  the  mech- 
anism of  the  various  interactions  and  the  ability  of  sug- 
gested models  to  fit  experimental  data.  Baronshkov  has 
produced  the  most  recent  data  on  mechanism  modeling  (Ref 
10).  Fig.  shows  a plot  of  the  average  nvunber  of  cascade 
particles  per  beam  proton  interaction  versus  the  energy 
of  the  incident  proton. 

For  copper,  there  are  approximately  nine  particles  pro- 
duced for  each  incident  5 UeV  proton.  The  qualitative 
break-out  of  these  nine  particles  is  as  follows: 

TABLE  I Distribution  of  cascade  particles  (?GeV  Beam) 

Yield  (E)  T/TOGeV  Y(t)  GeV 


Protons  (S) 

.7 

.39 

.27 

Protons  (C) 

2.0 

.015 

.03 

Neutrons  (S) 

• 7 

.39 

.27 

Neutrons  (c) 

2.0 

.015 

.03 

1.3 

.09 

.12 

n- 

1.3 

.09 

.12 

n° 

1.3 

.09 

.12 

9.3 

.96  GeV 

i 


The  (S)  above  stands  for  shower  particles  and  the  (C) 
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The  number  of  cascade  protons  per 
incident  primary  proton  as  a function 
of  atomic  weight  of  target. 


Fig. 5 Average  number  of  cascade  particles 
per  incident  proton. 
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for  cascade.  In  the  cascading  process,  the  initial 
10  GeV  proton  may  undergo  a collision,  producing  a proton 
of  energy  only  slightly  less  than  the  incident  particle 
energy;  these  protons  are  still  capable  of  proving  further 
cascading  and  thus  are  termed  cascade  protons.  If  these 
cascade  protons  (C)  interact,  they  produce  secondary  or 
tertiary  protons  with  somewhat  less  energy  (see  Tables 
I and  II),  and  these  lower  energy  protons  are  called  shower 
protons  (S). 


For  energies  on  the  order  of  10  GeV  the  qualitative 
breakout  is  as  follows: 

Table  II  Distribution  of  cascade  particles  (10  GeV  Beam) 


(E)  GeV 

Relative  # 

Energy  Fraction 

Primary  Protons  (C) 

8.5 

2x10" ^ 

1.7x10-^ 

Secondary  Protons  (S) 

2.9 

8.8x10" ^ 

2.6x10-^ 

Tertiary  Protons  (S) 

'\^1.0 

1.1x10" 5 

1.1x10-^ 

Secondary  Neutrons  (S) 

3.6 

8. 8x10" 3 

3.2x10-^ 

Tertiary  Neutrons  (S) 

'^1,2 

1.1x10" 2 

1.3x10"^ 

n+ 

.9 

2x10"'* 

1.8x10"® 

y+ 

.1 

6x10^ » 

6x10" ® 

As  can  be  seen,  the  mechanism  for  interactions  in  these  high 


energy  regions  are  continually  changing,  as  are  the  products. 

Fig.  5>  a plot  similar  to  Fig.  4,  extends  the  energy  of 
the  incoming  particle  to  show  a maximum  at  approximately  5 
to  10  GeV.  This  deviation  from  the  expected  (the  higher 
the  energy  of  the  incident  particle,  the  more  cascade  par- 
ticles) is  attributed  to  a change  in  the  density  of  the  nu- 
clear nucleons  during  the  course  of  development  in  the  inr 
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!i 

ternuclear  cascade  (Ref  9)«  The  large  number  of  shower  par- 
ticles produced  in  the  first  inelastic  N-l^or  N-n  collisions 

ij  inside  the  nucleus  knock  out  the  next  nucleons  along  the 

\ ! 

I path  of  cascade  development,  and  because  of  the  relativistic 

speeds  under  consideration,  the  nuclear  density  does  not 
I have  time  to  become  redistributed  appreciably  during  the 

i time  of  passage  of  the  nuclear  cascade,  or  simply  the 

nucleus  does  not  have  time  to  redistribute  its  nucleons  to 
! a state  which  would  allow  further  interactions.  The  cas- 

cade particles  encotmter  smaller  numbers  of  internuclear 
nucleons j thus,  the  excitation  energy  is  accordingly  smaller, 
and  consequently,  the  number  of  evaporated  particles  (de- 
layed emission)  is  lower.  The  energy  at  which  this  pheno- 
! : menon  occurs  is  termed  "the  saturation  energy"  and  may  pre- 

i 1 

; sent  a limit  on  the  energy  of  the  protons . 

I Let  us  now  look  at  some  of  the  secondary  particles  pro- 

t duced  in  the  cascade  process  in  order  to  determine  if  rea- 


sonable amounts  of  energy  may  be  deposited.  The  shower 
protons  have  somewhat  less  energy  than  the  cascade  protons, 
but  their  energies  are  still  high  enough  that  the  cross 
sections  are  low,  as  can  be  seen  by  comparing  Table  II 
with  Fig  2.  The  cascade  particles  have  similar  energies 
to  the  Incident  protons  and  thus  also  very  low  cross  secticns. 
Energy  may  be  deposited  due  to  multiple  scatters  of  the 
above  particles.  With  each  scatter,  the  cross  section  for 
the  next  collision  becomes  larger  until  the  particle  dis- 
sipates large  amounts  of  energy  with  each  collision. 
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Neutrons,  being  similar  to  protons  in  mass  and  energy- 
received  during  cascading,  have  very  low  cross  sections  for 
inelastic  interactions.  If  fissionable  material  were  present 
in  the  target,  the  increase  in  the  fission  barrier,  due  to 
deeper  splitting  at  higher  energies,  could  not  be  compen- 
sated for  by  an  Increase  of  excitation  energy,  and  the  fis- 
sion cross  section  should  decrease  with  increasing  energy  of 
the  primary  particle.  In  other  words,  if  a nucleus  is  to 
fission,  it  must  overcome  a barrier  caused  by  an  increase  in 
surface  tension  energy  of  the  two  new  nuclei  being  formed. 
This  increase  forms  an  energy  barrier  which  must  be  overcome 
by  the  energy  of  the  incident  particle.  The  lowest  poten- 
tial barrier  occurs  when  nuclei  of  the  same  mass  number  are 
formed,  but  with  higher  energy  incident  particles,  the  new 
nuclei  may  not  be  symmetric  and  thus  the  fission  barrier  is 
much  higher  with  this  deeper  splitting.  This  effect  has 
been  observed  for  plutonium  and  uranium  nuclei.  Neutrons, 
like  protons,  engage  in  nonfissionable  nuclear  interactions 
causing  partial  dissipation  of  their  energy  and  thus  continu- 
ing the  cascade  process. 

Large  number  of  plons  ^TT  are  produced  in  the  p-p 
and  p-n  reactions  of  the  internuclear  cascade.  These  par- 
ticles have  a mass  of  '\>135-1^0  MeV  and  possess  very  short 

half-lives.  The  neutral  pion,  n°,  has  a half-life  of 
16 

1.8x10  sec.  decaying  to  two  annihilation  gamma  rays. 

These  gammas  will  be  quickly  stopped  (a  drop  in  intensity 
of  '/e  in  2 cm  Al)  in  most  heavy  materials  due  to  the  high 


absorption  coefficient.  The  mechanism  for  absorption  will 

be  almost  totally  pair  production  and  show  absorption  co- 

-1 

efficients  near  .8cm  or  a drop  in  intensity 

where  Xo  - initial  intensity  of  the  radiation  avd  X - the 

intensity  after  a distance  X. 

For  a 70  MeV  gamma  I/Io  would  be  10  after  traveling 

through  one  meter  of  alumlniim.  The  7T±  have  considerably 

_ s 

longer  half-lives,  2.85  10  sec.,  and  since  they  are  mov- 
ing at  the  speed  of  light,  may  not  deposit  their  energy 
in  the  target  material.  The  rr+  may  decay  producing  a muon 
and  a neutrino.  The  muon  may  then  Interact  and  also  decay. 

77- j-  ^ 

e*  *iy  ■»  V' 

where  neutrino  and  v- = antineutrino. 

The  may  also  undergo  reaction  with  the  nucleons  and 

produce  a multiple  of  subnuclear  particles.  Of  the  energy 
3 

give  to  the  7T  , fairly  large  amovints  will  be  deposited  in 
the  target;  thus,  the  nvimber  of  pions  produced  becomes  a 
major  factor  in  determining  the  energy  deposition  in  a tar- 
get. 

Fig.  6 shows  the  average  number  of  pions  emitted  from 
the  various  elements  when  bombarded  by  protons  of  various 
energies.  Fig.  7 shows  the  distribution  of  the  pions 
between  77**,  77"^,  and  TT' at  various  energies  in  copper. 

The  TT^becomes  the  dominant  particle  produced  at  high  ener- 
gies and,  coupled  with  its  short  half-life,  would  Indicate 
a slight  increase  in  deposition  energy  from  this  process. 
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Nuclear  interactions  account  for  a major  portion  of 
the  energy  deposited  in  the  target  material.  However,  be- 

I 

cause  of  the  complexity  of  the  interaction  mechanisms, 
large  computer  programs  using  Monte  Carlo  techniques  are 

i 

used  at  the  present  time.  Computer  codes  such  as  HETC 
(High  Energy  Transport  Code)  (Ref  7)  do  not  take  into 
account  secondary  and  tertiary  particle  interactions.  Bar- 

I 

enshkov  has  recently  developed  a code  which  does  take  into 
accoiint  the  secondary  interactions  of  the  cascade  and  shows 
fair  agreement  with  experimental  resixlts.  Considerable 
modeling  and  experimental  work  must  still  be  completed  in 
this  area  before  application  of  theory  to  experimental  data 
involving  nuclear  cascading  will  be  completely  understood. 
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Bremsstrahlimg 

When  a proton  interacts  with  the  coulomb  field  of  the 
nucleus  of  the  target  material  through  ineleastic  collisions, 
the  accompanying  deceleration  leads  to  the  emission  of 
electromagnetic  radiation  whose  amplitude  is  proportional  to 
the  change  in  velocity  or  direction.  Bremsstralung  pro- 
duced by  collisions  of  fast  particles  with  the  atoms  will 
make  up  a continuous  spectriun  of  discrete  photons  whose 
wavelengths  are  related  to  the  energy  loss  of  the  particle 
(Ref  3): 

If  --  k-k' 

where  K and  K.'  are  the  energies  of  the  particle  before  and 
after  interaction  with  the  coulomb  field,  and  h equals 
Planck's  Constant,  C equals  the  speed  of  light,  and  ^ 
equals  the  wavelength  of  radiation  emitted. 

For  a small  deflection  in  a point  coulomb  field  of 
charge,  5 momentum  change  on  a particle  of  charge 

and  mass  (M)  is  transverse  and  is  given,  as  was  shown 
in  the  ionization  case  by: 

AP  - £ r f f — (/3) 

where  b is  the  impact  parameter  or  distance  of  closest  ap- 
proach to  the  nucleus  (Ref  3)* 

For  a relativistic  case,  if  we  transfer  reference 
frames  from  the  laboratory  frame  where  the  particle  moves 
to  the  frame  where  the  particle  is  at  rest,  classical 
theory  remains  accurate  with  only  a few  minor  changes  in 
the  impact  parameter.  Thus,  from  EQ  (I3)  and 


2k 


r 


V = ^ 

b 


The  intensity  of  the  radiation  can  be  shown  to  be  pro- 
portional to  the  square  of  the  product  of  the  amplitude  and 
the  charge  will  vary  as: 

r <i  ^ ^ i-  - 

/V\* 

Thus,  the  total  bremsstrahliing  varies  as  the  square  of  the 
atomic  number  of  the  stopping  material  and  Inversely  with 
the  square  of  the  mass  of  the  incident  particle.  Therefore, 
a proton  will  produce  about  one-millionth  the  bremsstrah- 
lung  of  an  electron  of  the  same  velocity.  Because  of  this 
strong  dependence  on  mass,  bremsstrahlung  is  almost  com- 
pletely negligible  for  all  fast  particles  other  than 
electrons . 
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Summary  of  Colllsional  Interactions 

In  summary,  let  us  look  at  examples  of  a few  of  the 
major  types  of  interactions  for  specific  cases  and  estimate 
the  energy  deposited  in  the  target  material.  These  cal- 
CTilations  are  only  to  the  zero  order  and  are  intended  to 
show  trends  and  not  meant  to  be  exact. 

As  a trial  case,  let  us  use  the  following  criteria: 
The  Beam 

1.  Energy  of  protons  = 10  GeV 

2.  (Xirrent  = 1 x Itf  amps 
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3.  Velocity  of  beam  = 2.98  x 10®“/s 
Relativistic  constant  = .988 

5.  Radius  of  beam  = .001  m 

28 

6.  Flxix  of  protons  in  beam  = 1.99  x 10 
protons/M^  sec 

_ 6 

7.  Charge  per  nanosecond  pulse  = 1.0  x 10 
coul/pulse 

The  Target 

2 7 

1.  Aluminum,  1 3AI 

2.  Mean  excitation  potential,  I = 166eV 

3.  Density  correction  faction 

§ = ^■.606X  + C = 4.02 

4.  Cross  section  for  interaction  = 11^ 
millibarns  = 115  x 10  m^ 

5.  Thickness  = 1 m 

The  calculations  will  be  done  for  each  mechanism  as 
follows : 

1.  Ionization  losses  will  be  calcxilated  via  Eq  (12) 
where  it  has  been  assumed  that  the  beam  does  not  slow  down 
with  energy  loss.  If  the  target  is  thick  with  respect  to 
the  range  of  the  proton,  this  is  not  a bad  approximation. 
For  this  example,  the  range  of  a 10  GeV  proton  is  on  the 
order  of  a few  meters;  therefore,  the  amount  calculated 

Is  somewhat  below  the  actual  value. 

2.  Nuclear  interaction  losses  are  calcvilated  on  the 
assumption  that  each  interaction  deposits  its  total  energy 
locally.  This  is  not  the  actxial  case,  since  we  discussed 


I 
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showering  in  a previous  section.  This  approximation  does 
not  take  into  accoiint  the  nuclear  energy  liberated  in  the 
interaction  processes  and  will  therefore  be  somewhat  low. 

3.  Bremsstrahlung  losses  will  be  estimated  from  the 
equation 

cL£.  E A^e,  / 

d E.  <»ll  ^ 7T  /37  ^fAr-ricl<  \C/  3 

which  relates  the  energy  loss  due  to  radiation  of  brems- 

strahl\xng  to  the  energy  loss  due  to  colliS('o/V5  (Ref  2). 

Table  III  clearly  shows  that  the  amoiant  of  energy 

deposited  by  the  beam  in  a target  is  surprisingly  small 

16 

compared  to  the  sublimation  energy  of  aluminum  (7*^xl0 
MeV/gram) ; thus,  there  is  some  concern  as  to  its  useful- 
ness. For  the  electron  beam,  we  would  expect  the  total  en- 
ergy of  the  beam  to  be  deposited  within  a very  short  dis- 
tance into  the  target  by  collisional  interactions. 

Table  III 

Summary  of  Energy  Deposited  by  Various  Interactions 


Mechanism  Max  Energy  Deposited  $ of  Total  ' ^ 

Ionization  3*5  x 10®  MeV/sec 

Nuclear 

Interactions  22.0  x 10®  MeV/sec 

Bremsstrahlung  0.6  MeV/sec 
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III.  Collective  Interactions 


Because  the  plasma  consists  of  extremely  large  numbers 
of  particles,  It  becomes  Impossible  to  follow  each  Individ- 
ual particle  and  Its  Interactions.  A more  practical  approach 
Is  to  represent  the  plasma  components  through  a distribution 
function.  If  the  kinetic  equations  are  then  applied  to 
this  function,  Instead  of  Individual  particles,  a "mean" 
response  will  evolve.  This  mean  response  can  be  considered 
the  response  of  the  entire  plasma,  or  a collective  response, 
and  is  felt  over  large  distances  as  compared  to  collisional 
effects . 

Oscillations  are  probably  the  most  basic  collective 
response  of  a plasma.  If  a plasma  Is  acted  upon  by  an  ex- 
ternal force,  it  will  tend  to  oscillate.  These  oscillations 
may  be  either  stable  or  unstable  In  nature.  If  the  oscil- 
lations are  unstable,  energy  from  the  force  will  be  con- 
tinually pulled  into  the  plasma  to  support  this  instability, 
and  this  energy  will  be  converted  into  random  motion  or 
heat.  This  is  the  mechanism  by  which  energy  from  the  pro- 
ton beam  can  be  transferred  Into  the  plasma. 

There  are  two  major  types  of  instabilities,  velocity- 
space  instabilities  in  which  the  plasma  is  unstable  be- 
cause of  its  particle  distribution  in  phase  space,  and 
configuration-space  instabilities  in  which  the  plasma  is 
unstable  because  of  its  particle  distribution  in  config- 
uration (or  real)  space. 

The  two-stream  instability  is  a well-known  instability 
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associated  with  charge  particle  beams  and  will  be  dis- 
cussed as  an  example  of  the  method  of  deposition  of 
energy  in  a plasma  through  unstable  oscillations. 

Two-stream  Instability 

Two-stream  instability  is  an  electrostatic  phenomenon 
involving  the  coupling  of  oscillations  within  the  electric 
field  of  the  beam  to  oscillatory  modes  within  the  plasma. 
Because  this  coupling  is  an  electrostatic  phenomenon,  and 
involves  separation  of  charge,  its  magnitude  is  much  larger 
than  other  electromagnetic  effects.  From  the  first  two 
moments  of  the  Vlasov  equation  for  the  beam  and  the  target 
plasma,  we  obtain: 

where  <^  = 1 for  the  beam  and  = 2 for  the  plasma  and 

f/i  (ii*) 

where  a = an  acceleration  term  due  to  forces  within  the  beam 
and  plasma 

N = number  density 
U = velocity  of  charges 
= energy  term. 

t 

Coupling  the  above  two  equations  and  solving  for  a yields : 

By  linearizing  the  above  equations,  they  become; 


Na  - Nt  -t- 

/V*' 

(the  1 for  superscript 

Udi,  s ^ 

means  a small  pertur- 

bat ion  term.) 
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thus  Eq  (15)  becomes 


(a(T  A/;  )(u^  + u;^) 


and  Eq  (16)  becomes; 

+ Nl){ut  4 (a/*  -+■ 

n:)(-hw-^K')  + (.Nt-^o) . (d) 

Keeping  only  the  first  order  terms,  Eq  (18)  becomes: 

^ ■+ 

Eq  (19)  becomes: 

4 Aiv/)  - j^{zf^iuiuu  X 

where  qi  = +e  and  q2  = -e,  and  there  is  no  zero  order  elec- 


tric field 


//  N?  4 < = O. 


With  the  above  equations,  we  have  h equations  and  h \inknowns, 

and  the  problem  can  be  solved.  If  we  wish  to  look  at  the 

normal  modes  of  the  two  systems,  assume; 

/w;  = A<,' 

where  and  u^  are  constants  Eq  (15)  then  becomes: 

^ - c K CN^ul  + N<uX)  . 

Eq  (16)  becomes:  (after  algebraic  simplification): 

= -<K(u:uJi)-h  (A//-/V';). 

Solving  the  above  four  equations  yields: 


. ) r o u;U*r.«,  ~ , fZl) 

^ /V\£o 


In  order  to  have  a nontrivial  solution  (u  0),  the  quan- 
tity in  the  brackets  must  vanish.  Therefore, 
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' (uJ-KuSy 

This  equation  will  have  four  roots  to  its  solution.  For 
example! 


Fig, 8 Solution  to  two-stream  instability  model (real  roots). 
Sketch  of  as  a function  of  wave  number (k). 


IThe  circles  represent  four  real  solutions;  however,  as 
gets  close  to  » in  other  words,  the  oscillations  of  the 

beam  and  plasma  couple,  we  get  two  real  and  two  imaginary 
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Two  a.OTJ 

Tw, 


k=u/u®  k=£i)/u® 


Fig,  9 Solutions  to  two-stream  instability  model  (imaginary 


Analytic  Estimate  of  Energy  Deposited  Through  Collective 
Interaction  in  a 1-D  Model 

As  can  be  seen  from  the  conclusions  in  the  section  on 
interactions , collisional  losses  do  not  deposit  large 
amo\ints  of  energy  in  the  target  material;  therefore,  the 
analysis  of  collisionless  effects  becomes  important  if  it 
can  be  shown  that  reasonable  energy  can  be  deposited  in  the 
target  by  these  effects.  As  a simple  test  of  the  procedure, 
the  following  first-order  analytic  solution  was  evaluated. 

If  we  linearize  the  first  two  moments  of  the  Vlasov 
equation  in  one  dimension 

J.  P _ __  A {/°U) 

■£  o(  X 

and 

--  - +/■“ 

where  p - particle  density 

u = velocity  of  particles 

P = pressure  term 

a = force  term  = § E 

m 

and  the  adiabatic  gas  law  which  is  a special  case  of  the 
third  moment, 

PoP 

where  K = i ^ s/f  / = the  degrees  of 

freedom.  By  letting 

P - P,  a--  a,4  a, 

u.  - u,  * U,  p ^ Po  -p, 

where  the  zero  order  terms  are  constants  and  the  first 
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order  terms  are  small,  the  resulting  equations  are 


r 


AA 

i-t 


P 


fi  ^ ' 
X 


^ CL 

J<X  ' 


(^3  ') 


^2*J) 


where  ol  - electron  speed  of  soimd,  U.  - o > and  all 

second  order  terms  are  neglected. 

The  result  is  Eqs  (zs),  (2*^)  and  (zs")  are  three  equations 
in  the  three  unknowns  p,u,*n^~P  • Solving  for  u,  take  the 
second  derivative  of  Eq  (zh)  with  respect  to  time 


/)  - 
^ <i-L 

- /f 

c^a, 

but  ^ = zz*"/?  so 

that 

/)  sllK 
^ Jt 

_ _ »(  a’- 

t. 

fi  Ja>  . 
<A-t 

\ 

1 

diXU-i 

Substituting  in  Eq  (23) 

- _ 


p 

and  canceling  the  yields 

. ^1-  JV, 


^ /2  ^ 


z a 


J<  X" 


ota. 

<^-t 


(2t) 


since 
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= the  field  due  to  the  beam 


where 

and  = the  field  due  to  the  plasma  motion. 

<La>  - 

«(t  f . y " e.  J 


since  jr  gu.fi,  and 

ified  value  obtained  from  outside  methods. 

A£±  - A u -i-  . 

A-t  ' IE 

Combining  Eqs  (2t)  and  (sn)  yields 


will  be  a spec- 
Thus , 

(2  7) 


a,  = _£.  iEl"  081 
7^.  ^ 


A -i  JL  )T  a-t 

Let  us  call  f A ^ , a forcing  ftinction,  and  note  that 

ic 

f >v»  = ^ ’■ 


/ /^  - / A/. 

/'<'l  f. 


A1  £■« 


where  Wp  is  the  plasma  frequency.  This  frequency  is  the 
frequency  at  which  a plasma  oscillates.  If  an  electron  is 
displaced  from  the  ion  site,  it  will  experience  a force 
that  seeks  to  return  it  to  the  original  position.  When 
it  arrives  at  its  equllibriiim  position,  it  will  possess 
kinetic  energy  equal  to  the  potential  energy  of  the  initial 
displacement  and  will  go  past  the  equilibrium  position,  thus 
setting  up  an  oscillation  of  frequency  mp.  Equation  (ie) 
then  becomes 


(2<)) 


Eq  (2*1)  is  a linear  differential  equation  in  Ui  for 
which  we  would  like  to  find  a solution.  Using  Green's 
functions,  if  F has  the  form 


F --  i(x-x')  I (i-f) 


then  the  solution  to  the  equation  is  Green's  function,  G 

u,  - 't') , (30) 

Then,  G must  satisfy  the  equation 

jtlL  -cC-dl^  + u// B = $ CX-x')  Sc£-dO.  (31) 

dl  d*- 

We  must  now  find  the  ftinctlonal  form  of  G.  cX,  for  any  F is 

^ (X-  X',  ^-t‘)  F (x\±‘)  c/F<U  i' . Csz) 

Multiplying  Eq  (ss'*  by  (~j^  ' +^’’)  yields 

^ i ^ ^ 

but,  from  Eq  (3/)  F (F,  F)  e(x  Ut' 

^ S a x-r')  S(±-^'>  FCxlF)  UF J-t ' 
then  for  any  G which  satisfied 

G ^ -d'}  = liX'X')  I H-t') 

then 


U,  - F(X,  t ) 

Thus,  we  need  to  find  a G for  which  this  fact  holds. 
Let's  take  the  fourier  transform  of  ^ (X^-X')  S(t-t') 

S(y-F)Ui-e)  - i^\oiKci.  (») 


let 


Gcx-x't-t')  = . 

(3^) 
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If  we  operate  on  Eq  (3‘^)  with  ) 


we  obtain 


JJ  A*-K*-+u;/)c 

(3S) 

We  would  like  Eq  (33)  and  (35")  to  be  equivalent,  thus 


where  0;*^=  is  known  as  the  Bohm-Gross 

dispersion  relation  and  relates  the  frequency  of  the 
electrostatic  wave  to  the  vector. 

A dispersion  relation  describes  the  ability  of  a spec- 
ific pulse  shape  to  pass  through  a plasma  without  being 
distorted.  A pulse  which  is  made  up  of  several  frequencies 
will  be  dispersed  upon  passing  through  a plasma  since  each 
frequency  travels  at  a different  speed.  The  extent  of 
this  dispersion  is  a fxmctlon  of  the  plasma  and  the  number 
and  range  of  frequencies  in  the  pulse.  Now  that  we  know 
j we  can  work  back  to  find  [j^x-x',  t-t') 

But  an  easier  method  may  be  to  do  the  following;  Let's 
look  at  the  fourier  transform  of  k (x.t) 

The  inverse  transform  is  then 

eft  (37) 


but  u(x,t) 


is  given  by  Eq  (32 ) thus 


[jCX-XU-^')  r(xU')cr""‘€'‘^^Ux'c(i:'  cCk  cL-L  . (38) 


Transforming  F( x\  t') 


and  Gfx-Xjt-t'J  yields 
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• * ' . / ^ ' 

^(xW}  = ™ /r^A':  W e e‘‘ 


' /\  I r r J ''  / “ f " »\  ' ^ ^ .€  a ~ ^ J 

tr  x-x,d-tt'>  ‘^ye  € . 

(tjo) 


Substituting  Eqs  (^o.)  and  (39)  into  Eq  (38)  yields 

u(ir.>^)  - rsiifni /T^v  cIk'ui' oCk' o(<^^  uk' ci<^' 


_ „ . ix'(K-kV 

f ( tf'.'-J')  L ( . •*^'')  C € 

e e 


J c/x 


<x(K-  k") 


X rK-K") 


etc.  so  that  where 


it  t ^ ft 


k'=  H"  , K”^  K,  UJ  : ui 


U(K.^)  --  F(K.tj)  3 , 


If  we  know  what  FC  X.w)  is,  then  we  can  find  Mfx.uj)  since 
we  also  know  G c «.  •*>) 
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If  we  model  the  beam  as  a chopped  system  in  which  the 

_i  1 

pulses  are  so  short  (10  seconds)  that  the  radius  of  the 
beam  and  its  length  are  approximately  equal,  the  pulse 
could  be  asstuned  to  be  a spherically  symmetric  charge 
distribution.  For  such  a sphere  of  radius  R,  the  electric 
field  at  a point  r is  given  by  (Ref.  3) 

F - ' 

when  r is  less  than  or  equal  to  R and  very  approximately 


equal  to 


for  r>R. 


c - _i_  J- 

^ ■ v7r£. 


Fig.  18  is  a plot  of  the  expected  force  (F  = qE) 
associated  with  the  model  of  the  pulse. 


Fig.  18  Force  due  to  spherical  charge  distribution  on 
a point  charge. 
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For  relativistic  proton  pulses,  the  electrons  first  see 
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i' 

i! 

; I 


I 

! 

I 


i 

I 


\ 

I 


1 
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the  effects  of  the  beam  at  fairly  short  ranges.  These 
forces  continue  to  grow  until  the  beam  just  reaches  the 
position  of  interest.  As  the  beam  passes  the  force  de- 
creases to  zero  at  the  center  of  the  pulse.  The  force  will 
then  Increase  to  a maximvim  at  the  opposite  edge  of  the 
pulse  and  decay  back  to  zero  as  the  pulse  moves  rapidly 
away. 


Figure  18  is  somewhat  similar  to  the  plot  of  the  de- 
rivative of  the  delta  function.  Thus,  the  force  on  the 
particle  was  taken  to  be  proportioned  to 

^ j ^ X-  v£  ^ 

where  c is  a constant  of  dimension  . 


A plot  of  ^ (x-v-6  ) 


would  look  like 


DirecT'OAJ 


where  first  the  electrons  are  pulled  toward  the  beam  and 
then  pushed  as  the  beam  passes. 

Transforming  Eq  ("f  I ) yields 
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Let  us  now  define  M = (Va). 


Where  M is  the  Mach  number 


of  the  beam  in  the  electron  gas . 


_ fc 


Then, 


d K . 


The  value  of  the  integral  in  Eq  (45)  will  depend  upon 
the  value  of  M,  the  Mach  number  of  the  beam  in  the  plasma. 
Solutions  to  the  integral  can  be  found  through  the  use 
of  Cauchy's  Integral  theorem,  the  theorem  of  residues.  The 
theorem  of  residues  states  that  the  integral  of  f(k)  around 
a closed  path  enclosing  a single  pole  of  f(k),  a point 
where  the  function  goes  infinite,  is  given  by  2Jl±  times 
the  residue  at  the  pole. 

dz  = 2ir/  (CL,,  + S-,  ) 


where  (X-x  represents  the  residue  of  the  pole  at  a.. 

Thus,  the  value  of  the  Integral  in  Eq  (45")  is  given  by 
2111  times  the  value  of  the  residue  at  any  pole  lying  within 
the  area  of  integration.  Let  us  divide  the  solution  to 
u(x,t)into  two  cases;  Case  1 where  M<1  and  the  beam  moves 
slower  than  the  speed  of  sound  in  the  plasma j and  Case 
where  M>1. 

If  M<1,  then  the  function  will  have  poles  on  the  imagin- 
ary axis  and  will  be  exponential  in  character.  Looking 
at  the  denominator,  we  find  that  the  value  of  the  poles 
is  given  when 


hi 


9BBB!2^B9S!^SC 


( 

goes  to  Infinity.  Thus  the  pole  can  be  fo\md  by 


/ 


I 


and  the  poles  are  located  at 

1 ,, 

a.(j-/y\'-)  '*■ 

on  the  imaginary  axis . 


(Mfe) 


plane,  and  the  residue  there  is  equal  to 

< uJp 


thus 
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UU,0  r 


Zl\cC-  O-M'-) 


(X-  , 

t -z^*-)  ''''- 


b)  If  (X-Mcct)<o,  then  the  beam  has  hit  and  passed 
through  the  point  of  Interest.  This  would  close  in  the 
lower-half  plaine  of  Fig.  13.  Thus,  the  only  pole  exists 
at  the  point 


K - 


thus 


^CX.e)  z ZjLICj^  e 

2 TT  Mg  a [ 

[ Z C-eco^)  . 

f / - /v\  »•)  V->-  ^ ' 


/CM 

Me  A _ /v\^  ) 


(47) 


where  the  limits  of  integration  require  us  to  add  an  addi- 
tional minus  sign  when  following  the  counterclockwise  rules 
of  residue  integration. 


If  we  now  plot  u(x,t)  ^<1, 


we  see 


aCZ-'^V 


X geAi^ 

b><Tc.nofJ 


Fig.l4  Solution 
to  intergral 


As  can  be  seen  from  the  equations  and  the  Fig.  l4, 
the  amplitude  of  u(x,t)  depends  upon  the  value  of  M.  The 
region  over  which  the  beam  is  felt  is  also  a f\mction  of 
M.  For  very  small  M,  exponential  fall  off  is  quite 
slow,  and  the  beam  can  be  felt  at  long  distances,  which 


is  as  would  be  expected.  As  the  beam  increases  its  energy 
and  the  Mach  number  increases  toward  1,  the  area  affected 
by  the  beam  becomes  quite  small.  In  other  words,  if  M is 
close  to  1,  the  point  of  interest  does  not  feel  the  beam's 
approach  until  it  is  almost  upon  the  point,  and  at  that 
time,  it  produces  a large  shock-front  type  Infxuence  on 
the  point.  Once  the  beam  passes,  the  beam  effects  decay 
rapidly.  Thus,  the  delta  function  beam  leads  to  a rapid 
acceleration  of  particles  at  the  points  close  to  beam 
passage . 

In  this  case,  M<1,  as  y goes  to  + or  - infinity,  u->-o 
because  of  the  shielding  effects  of  the  neutralizing  back- 
ground ions . 

As  a side  note,  if  we  let  M go  to  zero,  the  shielding 
distance  is  where  a.  can  be  written  as 


thus 


ujf  ^ \ 17^  J ' 


the  Debye  length.  The  Debye  length  is  the  length  at  which 
electric  effects  can  be  felt  in  a plasma. 


Case  II  where  M>1 


r [ JL^ 

XlTMeO.  ) {AA~i)'X'--  u)p'-/<\.'- 


The  poles  will  be  on  the  real  axis,  and  we  shovild  also 
note  that  since  the  beam  is  moving  faster  than  the  speed  of 
sound  in  the  plasma,  no  prior  information  can  be  passed 
ahead  of  the  beam.  Since  we  have  ignored  the  magnetic  field 
effects  of  the  beam  and  the  plasma,  only  the  longitudinal 
wave  will  be  present,  and  its  velocity  will  approach  the 
speed  of  sound  in  the  plasma  instead  of  the  speed  of  light. 

Fig.  15  shows  the  relation  between  the  phase  velocity 
of  a wave  with  magnetic  effects  and  plasma  wave  without 
magnetic  effects. 


Fig.  15  Phase  velocity  as  function  of  normalized  frequency 
for  transverse  and  longitudinal  waves . 

At  low  normalized  frequencies,  the  phase  velocities  are 


unboamded;  however,  as  the  frequencies  increase,  the  phase 
velocity  of  the  electron  plasma  is  smaller  than  the  velocity 
of  the  beam,  and  no  electromagnetic  effects  will  be  seen 
before  the  wave  strikes.  This  leads  to  the  fact  that  until 
t = o,  u(Xjt)  must  integrate  to  zero,  and  therefore  no  poles 
can  be  located  in  the  upper-half  plane  of  the  figure. 

Thus,  for  Case  Il-a  where  (x >0 

- 0. 

The  remaining  case  ( Il-b)  in  which  (x-amt)<0  would 
be  described  by  the  lower  half-plane  of  Fig.  16 


Fig.  16  Location  of  residual  poles 


where  the  poles  have  been  displaced  from  the  axis  by  ie,  the 
denominator  of  Eq  (q  U)  then  becomes 

(M^-0  + 


F-- 


III!  1111  t ipiwm — 


Thus,  the  poles  are  located 


taking  the  limit  as  e-*-0  yields  K = - 
Thus,  from  the  residue  theorem 


urx^o-  - 


TTi] 


a y«- 


-l  u)f  X 


t-  /)  / - „ ) 

V ^ a <yn‘-0 


if7 
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Fig.  17  shows  a plot  of  Case  Il-b  M>1 


where 

An  times,  prior  to  the  beam  hitting  the  point  of 
interest,  nothing  happens.  After  passage,  we  now  have 
waves  traveling  at  V = Ma  of  the  beam,  without  dispersion 
at  a frequency 


IfJ  - \jip  /v\ 

N 

(.M*--  1)'^- 

M 

> 

Vf  , 

*L  - Jf^ 

^ - 
iK  ' 

TKt.  5rc.c;f 

W 

From  Fig.  I7,  we 

can 

also  see  that 

the 

amplitude  of  the 

j wave  is  proportional  to  ^/M;  thus,  for  maximum  interaction 

we  would  want  v^i^a  or 


* 
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RS 


a I V.  y yL/o'\r^ 

^ tA  t ~ jrAt*^ 


= 66  X lO^m/sec  which  is  far  below  relativistic 


speeds.  Also,  as  M becomes  large  oj  -►  cUp,  and  as  expected 


the  plasma  will  not  oscillate  below  the  plasma  frequency. 

The  maximiam  electric  field  associated  with  the  pulse 
is  given  by 


f 


^ne.  <?*• 


f.  o ^ 


where  q = the  total  charge  in  the  pulse,  IxlO”  coul),  and 
R = the  radius  of  the  pulse  (.001  meters). 

From  Fig.  17,  the  amplitude  or  maximum  velocity  of  the 
electron  was  given  by 


fl  -- 


fC 


Alt  V->. 

where  - '^A  = M-515 


.16 


= = 2.4  X lO"*^  sec. 

fm 


C = the  forcing  function  constant 


gc  _ J: 

■ AV  4.  t 


.(XT 


Therefore,  C S'  = £ 

C = E X ( a.  radius  for  the  S'  function) 


C = E:  R'  r (fo  ^ 

12  Q. 

thus,  /?  = 4.3  X 10  /s  which  is  faster  than  the  speed  of 
light  and  the  model  breaks  down  at  this  point.  If  we  con- 
tinue the  analysis,  we  find  the  kinetic  energy  of  each 
electron  is  given  by 


£7  = '4AtcR'=  5.2  X 10®MeV 
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This  amount  of  energy  is,  to  say  the  least,  not  a 
small  perturbation,  and  the  linear  analysis  is  clearly  in- 
correct; however,  it  does  indicate  that  energy  of  consid- 
erable amoxmts  may  be  transferred  by  this  process,  and 
further  investigations  are  warranted. 

If  we  scale  the  charge  in  the  pulse  by  a factor 
11 

of  1000  down  to  10”  coulombs,  the  velocity  of  the  plasma 
electrons  is  equal  to  4.2  x 10’“^/s  and  relates  to  an  energy 
of  approximately  5.2  KeV.  At  this  level,  the  energy  of 
the  actual  beam  (10~  coulombs)  would  be  considerable  if 
scaled  up  to  the  10 


coulomb  value . 


IV.  Method 


As  was  demonstrated  in  the  previous  section,  there  is 
a need  for  further  understanding  in  the  area  of  configur- 
ation-space instabilities.  Configuration-space  instabili- 
ties occur  when  the  potential  energy  of  a system  is  lowered 
through  the  redistribution  of  its  spatial  components.  Thus, 
to  reach  this  lower  energy  level,  the  members  of  the  plasma 
must  flow  to  a new  spatial  location.  The  movement  of  par- 
ticles in  a plasma  is  best  described  by  the  Boltzmann  equa- 
tion. It  shows  the  change  in  the  distribution  of  particles 
in  a plasma  subject  to  both  internal  and  external  electro- 
magnetic fields  and  collisional  interactions.  For  this 
thesis,  we  have  assumed  that  the  collective  interactions 
take  place  on  time  scales  much  shorter  than  the  binary 
collision  times;  thus,  the  collision  effects  on  the  particle 
distribution  can  be  ignored.  This  approximation  is  commonly 
called  the  Vlasov  approximation,  and  as  expected,  the  one- 
particle  Boltzmann  equation  becomes  the  Boltzmann-Vlasov 
equation 

4-  /^S-*-(vx5))  Vp'f  = o 

where 


5 = 
■t  = 
V = 


the  one  particle  distribution  function 
time 

vector  velocity  of  the  particle 
the  Lorentzian  force  acting  on  the  particles. 
The  force  is  composed  of  the  electric,  ft  , 
and  magnetic,  <g(vKi^  forces  from  both  in- 
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ternal  and  external  sources . 

The  development  of  the  Vlasov  equations  from  the 
basic  equations  associated  with  particle  motion  (Liou- 
ville's  equation)  involves  many  assumptions  that  are  not 
readily  apparent;  thus,  the  development  of  this  equation 
is  done  in  Appendix  A. 

Exact  solutions  to  the  Boltzmann- Vlasov  equations 
become  very  difficult,  because  the  distribution  functions 
are  defined  in  phase  space.  Large  amounts  of  computer  stor 
age  would  be  necessary  to  describe  any  distribution  in 
phase  space  and  would  become  practically  impossible  to 
follow  during  interactions . 

Thus,  the  moments  of  the  Vlasov  equations  were  taken, 
and  the  velocity  space  was  Integrated  out,  leaving  the  equa 
tion  in  real  space  and  a more  tractable  problem. 

Through  the  moments  , the  variation  in  density,  average 
momenttim  of  the  particles,  and  average  temperature  or  press 
ure  can  be  readily  computed.  Other  physical  parameters  of 
the  plasma  can  be  investigated  by  taking  higher  order 
moments;  however,  the  complexity  of  the  equations  increases 
rapidly  beyond  the  third  moment. 

The  first  three  moments  of  the  Vlasov  equation  are 
usually  termed  the  hydro - equations , because  they  involve 
the  flow  or  redistribution  of  the  particles  and  the  physi- 
cal properties  associated  with  that  flow.  Thus,  we  shall 
develop  the  first  three  moments  of  the  Vlasov  equation  to 
examine  the  response  produced  by  a charged  particle  beam 


in  a target  plasma. 


I 
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For  our  moment  analysis,  two  models  will  be  used:  a 

one -dimensional  model,  and  a three-dimensional  model  with 
two  degrees  of  freedom,  hereafter  to  be  referred  to  as 
two-dimensional . 

The  one -dimensional  model  was  used  to  establish  the 
procedures  for  modeling  the  target  and  the  beamj  and  served 
to  iron  out  major  problems  which  develop  in  computational 
and  physical  analysis.  The  beam  target  interaction  area 
closely  resembles  a cylindrical  volume;  thus,  the  two-di- 
mensional model  will  be  structured  in  cylindrical  coordi- 
nates with  angular  symmetry.  As  the  beam  strikes  the  tar- 
get, it  is  expected  that  there  will  be  variations  in  den- 
sity, momentum  and  temperature  in  the  axial  and  radial 
directions,  but  none  in  the  angular  direction.  Thus,  the 
assumption  of  angular  symmetry  reduces  the  number  of  de- 
grees of  freedom  to  two.  Other  specific  assumptions 
will  be  made  and  discussed  with  each  model,  but  the  pro- 
perties of  the  beam,  plasma  structure  and  solving  routine 
will  be  common  to  both  models  and  discussed  before  the 
specific  models  are  introduced. 

The  Plasma 

The  target  material  is  assumed  to  be  an  aluminum  slab 
of  approximately  one  meter  in  thickness.  The  outer  shell 
electrons  in  aluminum  (3S*,3PM  are  free  to  move  in  the 
conduction  bands  of  the  ionic  lattice.  The  assumption  is 
made  that  Fermi-Dirac  statistics  hold  for  these  electrons 
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or  that  they  are  indistinguishable  and  cannot  be  located  si- 
nultaneously  in  the  same  state.  As  the  beam  deposits  en- 
ergy in  the  plasma,  the  temperature  and  density  will  change. 

In  the  limit  as  the  temperature  increases  and  the  density 
decreases,  Fermi-Dlrac  statistics  become  Maxwell -Boltzmann 
statistics,  and  since  our  simple  analytic  example  Indicates 
considerable  heating,  the  assumption  was  made  that  the 
plasma  could  be  represented  by  a Maxwellian  distribution. 

The  electromagnetic  forces  which  act  upon  the  metallic 
structure  will  tend  to  displace  both  the  electrons  and  the 
residual  ions,  but  because  of  the  large  mass  difference  be- 
tween the  electrons  and  ions,  the  electrons  will  move  much 
farther  in  a given  time  and  with  much  greater  velocity. 

Compared  to  the  electrons,  the  ions  appear  stationary  over 
short  time  scales. 

It  will  be  convenient  to  represent  the  target  as  an 
electron  gas  with  a Maxwellian  distribution.  The  ions, 
although  they  do  not  move,  serve  as  neutralization  sites 
and  fiinction  as  stationary  positive  charge  sites  for  the  re- 
storing electric  field  developed,  if  the  electrons  shoxild 
be  displaced.  Thus,  the  target  may  be  categorized  as  a 
single  species,  Maxwellian  distributed,  charge  neutralized, 
electron  gas. 

The  Solving  Routine 

The  solving  routine  used  in  this  thesis  was  an  eight- 
point  ordinary  differential  equation-solving  routine  which 
used  the  Newton-Raphson  method  of  solution.  The  routine 

i 
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solves  a system  of  first  order  ordinary  differential  equa- 
tions at  a point  B,  given  initial  conditions  at  point  A, 
and  a subroutine  for  determining  the  values  of  the  partial 
derivatives . 

Because  the  equation- solving  routine  can  be  extremely 
complex  or  extremely  simple,  depending  upon  availability 
to  the  user,  no  specific  emphasis  will  be  placed  on  the 
solving  routine  in  this  report. 

The  Beam 

The  proton  beam  is  assumed  to  be  a 10®  amp.  beam  of 
10  GeV  protons  chopped  into  pvilses  at  nanosecond  intervals. 
Because  of  the  poor  propagational  properties  of  charged  par- 
ticle beams  in  the  atmosphere,  considerable  amounts  of  en- 
ergy per  particle  and  fairly  high  density  will  be  required 
to  propagate  the  beam  over  any  reasonable  distance.  Ac- 
celerators capable  of  producing  such  high  energies  and  den- 
sities are  not  yet  available;  thus,  these  figures  should 
be  considered  only  as  representative  of  actual  values.  With 
such  high  energies,  the  properties  of  the  beam  will  clearly 
have  to  be  treated  in  a relativistic  manner,  and  the  range 
of  these  particles  in  solids  will  be  on  the  order  of  a few 
meters.  As  was  discussed  in  the  collisional  interaction 
section,  the  energy  lost  by  the  beam  will  take  place  over 
the  entire  range  of  its  interaction,  with  ionization  losses 
becoming  important  only  near  the  end  of  the  particle's 
range. 
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We  have  assumed  the  target  to  be  one  meter  thick;  thus, 
let  us  also  assume  the  energy  lost  by  the  beam  in  passing 
through  the  target  is  small  and  that  the  beam  does  not  slow 
down.  Since  we  are  more  concerned  with  the  effects  on  the 
plasma  caused  by  the  fields  of  the  beam,  this  assumption 
should  not  induce  any  major  errors. 

The  chopping  of  the  beam  into  pulses  was  assumed  be- 
cause of  the  hose-instability  phenomenon  which  occurs  in 
beam  propagation.  The  instability  so-called  because  the 
beam  will  tend  to  act  like  a high-pressure  fire  hose  if  not 
restrained.  By  pulsing  the  beam,  a major  portion  of  this 
instability  can  be  eliminated.  If  the  pulses  are  short 
enough  and  we  assume  the  radius  of  the  beam  to  be  small  in 
relation  to  the  target,  one  could  consider  the  beam  to  be 
a series  of  point  charges  pulsed  at  nanosecond  intervals. 
This  has  been  the  assximption  made,  and  the  electric  field 
associated  with  the  model  is  calculated  in  Appendix  B.  The 
magnetic  field  of  the  beam  has  been  ignored  in  both  models. 
In  the  one-dimensional  model,  the  field  is  not  defined 
since  the  dimension  of  choice  is  that  along  the  direction 
of  travel  of  the  beam.  In  the  two-dimensional  model, 
the  magnetic  field  will  produce  an  outward  radial  force 
on  the  plasma;  however,  for  the  initial  look  into  col- 
lective interactions,  only  the  electric  field  will  be  con- 
sidered. Additional  future  work  will  have  to  include  the 
magnetic  forces  of  the  beam. 
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The  One-dimensional  Model 

In  the  development  of  the  models,  the  space  under  con- 
sideration will  first  be  discussed.  Because  of  the  large 
areas  under  consideration  and  the  small  areas  of  interaction, 
scaling  parameters  will  be  applied  to  the  various  parameters. 

The  development  of  the  moments  will  then  be  undertaken.  | 

The  internal  and  external  electric  and  magnetic  fields  will 
be  coupled  to  the  moments  and  initial  conditions  specified. 

One-dimensional  space.  The  one-dimensional  model  | 

assiunes  variation  only  in  the  direction  of  motion  of  the  | 

i 

beam,  and  since  the  solving  routine  uses  an  eight-point 
fitting  routine,  the  number  of  nodal  points  (points  of  solu- 
tion) will  be  multiples  of  eight,  plus  an  additional  node  ; 

to  accovint  for  the  Infinity  of  space  outside  the  target  I 

material  even  though  no  solution  will  be  obtained  at  this 
node.  The  beam  will  initially  be  at  X=-®  at  time  t = j 

and  will  move  toward  +oo  as  time  Increases.  Since  the  tar- 
get is  finite  in  dimensions,  most  of  the  time  the  beam  will 
be  too  far  away  to  be  of  any  consequence  and  of  little  in- 
terest. To  eliminate  this  problem  of  extremely  large  axial 
distances^a  scale  transformation  was  chosen  to  allow  move- 
ment from  X = - to  + infinity  over  the  scaled  range  of 
W = +1  to  -1.  By  picking  the  proper  function,  the  region 
of  interest  can  comprise  the  greatest  amount  of  scaled 
space,  and  the  areas  further  away  can  be  drastically  con- 
densed. 

The  function  Z - Z<,  ^ serves  such  a purpose. 
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The  constant  can  be  adjusted  to  concentrate  on 
areas  of  interest.  If  we  let  A u r '/<?  and  run  u between 


I 
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-1  and  +1, 


we  obtain  18  intervals . 


-%  -'h  -1^  o Ifif  ^ 34,  ^ 7/f  ^ 

I— I — t— L I I I I I I 1.  I I I I 1 I I I 

D,i-e<.rfo n/ 

If  we  let  zero  equal  the  center  of  the  target  and  the 
points  ± ®/9  equal  the  end  of  the  target,  we  will  concen- 
trate the  area  of  little  interest  in  the  last  cell  on  either 
end.  These  cells  will  not  be  analyzed  since  we  have  the 
17  inner  points  for  the  solving  routine.  This  method  also 
may  prevent  problems  with  fitting  the  large  jiunp  in  the 
density,  momentum  and  temperature  between  the  target  and 
no-target  regions. 

In  a similar  manner,  the  time  and  the  momentum  of  the 
electrons  were  scaled  to  concentrate  on  areas  of  interest. 

We  wish  to  examine  the  plasma  at  the  time  of  Impact  of  the 
beam  without  doing  much  computation  in  the  areas  of  little 
Interest.  The  scaling  function 

T = To  Sinh(s)  (SZ) 

was  selected.  This  fxinctlon  takes  very  small  steps  at  times 
near  zero  (when  beam  reaches  center  of  slab)  and  large  steps 
with  increasing  values  of  (s) . By  adjusting  To,  the  speci- 
fic times  of  interest  can  be  used. 

The  momentum  is  likewise  a candidate  for  scaling  with 

p=  siA/hCv)  (^-i) 
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Thus,  the  region  of  concentration  is  the  region  of 
small  momentum.  If  large  amounts  of  energy  are  deposited 
by  the  beam,  Po  can  be  adjusted  as  necessary. 

The  Moment  and  Force  Equations . In  one  dimension,  the 
Vlasov  equation  reduces  to 

-j£  (X,.)  * +5^=0.  fe.) 

The  force  term,  F = (vxa) , also  reduces  since  the  magnetic 
field  associated  with  the  beam,  and  any  internal  magnetic 
field  caused  by  charge  motions,  will  be  undefined  in  one 
dimentlon.  Thus,  the  Vlasov  equation  takes  on  the  form 


where  p = mv,  and  it  is  assumed  that  the  motion  of  the  elec- 
trons is  non-relativistic . If  considerable  energy  is 
dmaped  into  the  system  by  the  beam,  this  equation  will  re- 
quire the  momentum  to  take  on  the  form  p = if  V , and  a 
f^llly  relativistic  model  will  have  to  be  examined. 

The  external  electric  field,  E^x-t  is  developed  in 
Appendix  B,and  estimates  of  the  electric  field  due  to  the 
beam  at  specific  points  are  given  by  the  function 


Q r \ _ vt 

vTTf.  K7^J 

— - — — • 


(•S’  (ti) 


Where  Q = the  total  charge  per  pulse  in  the  beam  (See 
Appendix  C) 

= the  relativistic  constant 
V = the  velocity  of  the  beam 

-t  = the  time,  where  the  beam  is  at  the  center  of 
the  target  at  ^ = o 
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R = radius  of  interaction  - in  1-D  chosen  as  an  arbi- 
trary constant 

X = the  position  iznder  consideration  . 

The  internal  electric  field  is  obtained  from  the  solu- 
tion to  Maxwell's  equation: 

(yxB)  f.  jt  + j (s^> 

where  J is  the  current  density. 

Since  B=o  from  symmetry  and  the  current  density  is  given  by 

1 - ^ j P ^CP)  (^9  ■ 

Eq  (57)  simplified  to 

A^.  - -i—  { P i CP'i  A p . 

,kt  - e,/^  } ^ 

In  Appendix  D,  the  first  three  moments  to  the  one-di- 
mensional Vlasov  equation  and  the  equations  necessary  to 
solve  the  time  derivative  of  the  internal  electric  field 
are  developed. 

The  analysis  in  Appendix  D resulted  in  four  equations 
for  the  time  derivatives  of  the  density,  momentum,  temper- 
ature and  the  internal  electric  field. 

The  first  moment,  continuity  equation 

The  second  moment 


The  third  moment 


where  H = the  particle  density 

P = the  average  momentum  of  an  electron 

'n*' 

P = a modified  pressure  term  since 
Ad 

If  the  scale  factors  and  constants  are  now  applied  to  these 
equations,  they  reduce  to: 

The  first  moment 

-j^z  R L>i  +4^] 


The  second  moment 


The  third  moment 


- rL^x.#- + y.  ^4-J  (4  3) 

where  Yi  the  scaled  density  parameter 

/i.  = the  scaled  momentum  parameter 
Y3  = the  scaled  press\ore  parameter 
R = a series  of  condensed  constants 
E.  = an  electric  field  scale  factor 
^ = time  (scaled) 

^ = distance  (scaled)  , 

The  internal  electric  field  is  given  by 

r - 2 £.  Coshes')  y^.  cxP(/<)  . (^,‘^) 

^ 5 

Equipped  with  Eqs  (il),  (i2),  (^3),  and  (i4),  and  the  dif- 
ferential solving  routine  described  earlier,  the  value  can 
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! be  stepped  forward  in  time  once  given  a set  of  initial  con- 

I ditions . 

The  Initial  Conditions.  The  initial  conditions  for  the  one- 
dimensional model  are  as  follows;  j 

•I 

1)  The  plasma  is  initially  cold.  Thus  at  ] 

i t = 

it.  - o ^ the  momentum  equals  zero 
iz  = o ^ the  temperature  is  zero 

2)  ii  , the  initial  density  distribution  is  given  by 
a slightly  smoothed-out  step  function.  Fig  (10)  shows  the 
initial  distribution  of  particles  normalized  to  let  1 = the 
density  of  the  plasma  electrons. 

This  initial  distribution  was  obtained  by  use  of  the 
following  exponential 

= exP(-3Z  + iz  i.izs-u  + .oi-s’ns-') 

3)  At  time  ^ the  external  electric  field  should 

be  zero,  and  the  position  of  the  beam  pulse  should  be 

X = -CD  • 

J 

The  parameters  of  the  one -dimensional  model  are  now  | 

j 

specified.  The  one -dimensional  computer  program  can  be  ! 

found  in  Appendix  F.  ^ 

The  Two-dimensional  Model 

The  two-dimensional  model  has  many  of  the  same  features 
as  the  one-dimensional  model,  and  where  similar,  the  de-  I 

tailed  analysis  will  not  be  given. 

Two-dimensional  Space.  The  one -dimensional  model  was  2 
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Fig.  10 

INITIAL  DISTRIBUTION 


a rough  approximation  to  the  real  physical  setup  of  the 
problem.  The  two-dimensional  model  is,  however,  a more 
realistic  model  of  the  problem.  The  beam,  considered  as 
a point  charge,  will  Interact  symmetrically  in  the  angular 
direction  of  the  cylindrical  coordinate  system.  In  the 
axial  direction  this  will  not  be  the  case  vmless  the  charge 
pulse  is  exactly  at  E = o , the  center  of  the  target.  For 
this  reason,  the  target  is  considered  a cylinder  of  finite 
radius  and  length,  being  struck  along  the  axial  line  of 
symmetry  by  the  point  charge. 


The  symmetry  associated  with  motion  of  the  plasma  in 
the  angular  direction  leads  to  the  assumption  that  there 
are  no  variations  in  density  and  momentum  in  this  devia- 
tionj  thus,  in  the  analysis,  these  terms  will  be  set  to 
zero. 

As  was  done  in  the  one-dimensional  model,  the  distance, 
momentum  and  time  parameters  were  scaled  in  order  to  con- 
centrate the  analysis  in  areas  of  greatest  interest.  The 
same  scale  factors  will  be  used  in  the  two-dimensional 
model  with  the  addition  of  a scale  factor  in  the  radial 
direction,  (^4, 

In  a similar  manner  as  was  done  in  the  axial  direction 


eh 


I 


i 

1 

i 


1 

! 

i 


1 

of  the  1-D  problem,  w can  range  from  o to  1 in  scaled  values 
corresponding  to  o and  infinity  in  real  space.  The  interval 
AW  will  also  be  '/<>  j thus  by  a proper  choice  of  (?„  the  last 
cell,  although  ignored,  contains  all  of  the  space  outside 
of  the  cylinder  in  the  radial  direction.  The  model  space 
therefore  looks  like  a planar  grid  with  17  divisions  in 
the  axial  direction  and  9 spaces  in  the  radial  direction. 


The  Moment  and  Force  Equations.  The  Boltzmann- Vlasov 

equation  in  tensor  notation  is  given  by 

«/f  + p^A£  ^ JL  = o /47) 

IT  JLX<  ^ 

where  equals  the  summation  over  1,  2 and  3 and  corresponds 
to  the  coordinates  iE 

f = the  particle  distribution 

P*^  = the  momentum  in  the  d direction 
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^ - the  force,  for  our  puipose  the  Lorentzian  force 
f C£*^+ ( VX  8)*^  , in  the  eK  direction. 

The  first  three  moments  of  Eq  (^7)  were  taken  to  yield. 

(See  Appendix  E): 

The  first  moment 

The  second  moment 

(ill 

* s!i 

The  third  moment 

•.  (7») 

where  /f  = Number  density  of  electrons  or  number  of  elec- 
trons per  cubic  meter  per  \init  momentiim 
/7  = The  momentum  in  the  i direction  where  1 = 1,3 
or  1 = r,  z,  respectively 

/i  = The  temperature  term  in  the  1 direction  where  Ij 
can  form  any  combination  1,  2 or  3 corresponding 
to  r,  4>,z.  (rz  = zr  etc.  in  Eulerian  theory) 

^‘  = the  force  on  the  electrons  in  the  1 direction. 

And  the  subscript  e<^ ) means  derivative  over  all  components 


AJL  + AJL  + 

jr  diZ  * 


66 


in  order  to  close  the  set  of  equations;  ijk 
can  be  any  combination  of  r,  41  and  z. 


After  adding  the  scaling  factors  and  making  the  simpl- 
fication  = o),  the  three  moments  reduce  to; 

The  first  moment 


The  second  moment 


a"  --  3^-  ^ i4 

* i fir- 
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ii  - f-r/f  J J (,s-) 


ii“  r (e')‘ ^(/iT 


TT 


+ r/?V  + + 


r C^TC/i'j  t f/f'/zf'’ 
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where  the  off-diagonal  components  of  the  pressure  in  the 
(j)  direction  have  been  set  to  zero  because  of  symmetry.  The 
momentum  component  in  the  ({>  direction  is  by  the  same  argu- 
ment zero. 

The  only  term  yet  unspecified  is  the  force  term.  The 
Lorentzian  force,  as  in  the  one -dimensional  model,  is  com- 
posed of  the  forces  due  to  both  external  and  internal  fields. 
The  external  fields  have  been  limited  to  only  the  electric 
field  from  the  beam.  In  the  two-dimensional  case,  however, 
there  is  both  an  axial  and  a radial  component  to  the  field. 
Appendix  C derives  the  equations  used  to  evaluate  these 
components.  The  axial  component  is  given  by  the  same  equa- 
tion as  in  the  one-dimensional  model. 


Q 

VTT£ 


V-t 


E{X)  : 


(Vi  ^ i^^d-  P")  J 


(78) 


where  = the  total  charge  in  the  point  beam 

^ - the  relativistic  constant 

^ = the  velocity  of  the  beam 

t = the  time 

f?  = the  radial  position  of  the  point  of  Interest 
X = the  axial  position  of  the  point  of  interest 
except  that  the  value  of  R will  not  be  a constant  but  will 
vary  with  the  position  of  the  point  of  interest. 

The  radial  component  of  the  electric  field  is  given  by 
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Q RCI-P*; 

V TT  €e 

- f/z  - /f?Vi-PV 

' (7<t) 

The  Internal  fields  consist  of  both  electric  and  mag- 
netic components  and  are  described  by  Maxwell's  equations 

IB  ^ t 

X t 


and 


IE 

<k  t 


X B 


X 


In  Appendix  E,  these  two  equations  are  expanded  in 
cylindrical  coordinates  and  scaled  to  the  parameters  used 
in  the  model.  The  time  derivatives  of  these  two  fields 
then  become: 


K 


C 

I 


(8  o') 


IB*  _ 

^ t ,A.  r/  • 


(8  1) 


(8^) 


The  force  in  the  radial  direction  is  then  given  by 
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J-r  - -eCe'’- 

and  In  the  z direccion  by 

= -e  C‘'J  f8“' 

Equipped  with  Eqs  (7/),  (7  2),  (73),  (74),  (72r), 

(74,),  (77),  (80),  (80j  and  (82)  and  the  solving  routine 
described  in  an  earlier  section,  the  values  of  the  density, 
momentum  and  temperature  can  be  stepped  in  time  from  a set 
of  initial  conditions. 

Initial  Conditions . The  initial  conditions  for  the 
two-dimensional  model  are  as  follows: 

1.  The  plasma  is  initially  cold.  Thus,  at 

= o for  i = r,z 

- e'T 

ft.  - o for  ij  = rr  rz,  (})()>,  zz 

2.  All  other  parameters  remain  the  same  as  the  one- 
dimensional  model. 
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V.  Results  and  Conclusions 


In  this  section,  the  sample  runs  which  would  estab- 
lish credibility  to  the  model  and  methods  used  will  be 
discussed  first,  followed  by  an  analysis  of  data  received 
from  both  the  one-  and  two-dimensional  models.  A short 
conclusion  section  and  a few  words  about  "what  should  be 
done  next"  will  follow. 

Validity  Runs 

The  most  basic  configuration-space  property  of  a 
plasma  is  a simple  oscillation.  The  plasma  system,  as 
described  by  the  first  three  moments,  resembles  a simple 
harmonic  oscillator.  V/hen  the  particles  are  initially  given 
a small  amount  of  momentum  in  one  direction,  the  internal 
electric  restoring  force,  which  develops  because  of  this 
displacement,  will  accelerate  the  particle  in  the  opposite 
direction  of  the  displacement.  As  the  particle  returns 
to  the  neutral  position,  it  possesses  some  kinetic  energy 
and  continues,  thus  the  process  can  be  represented  as  the 
perfect  spring  model  for  simple  harmonic  motion.  Appen- 
dix G is  a short  derivation  of  the  frequency  of  this  har- 
monic motion  and  shows  the  frequency  to  be 
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where  uJp  = the  plasma  frequency 

N = the  number  of  electrons  per  cubic  meter 
/'''«=  the  mass  of  the  electron 
€o  = the  permittivity  of  free-space 
If  our  model  could  demonstrate  these  oscillations  at 
the  proper  frequency,  then  some  credibility  must  be  given 
to  the  equations  and  techniques  used.  Table  IV  summarizes 
the  results  of  the  introduction  of  a small  amount  of  mo- 
mentum into  both  models. 


Table  IV  Summary  of  plasma  oscillations. 


Model 

(sec~^) 

Calculated  (Infinite) 

2,k  X 10 

1-D 

+ 15 

3.12  X 10^ 

2-D  axial 

+ 15 

3.8  X 10 

2-D  radial 

+ 15 

3.8  X 10 

Although  the  calculated  frequencies  listed  are  not  as  close 


to  the  calculated  values  as  one  might  like,  it  should  be 
noted  that  the  calculated  value  is  for  an  infinite  model 
and  our  model  is  finite.  This  fact  should  decrease  the 
value  of  (J<>as  will  be  shown  in  Appendix  G and  thus  make 
our  calculated  values  much  closer.  Additional  runs  were 
made  to  verify  the  analysis  by  varying  both  the  density  of 
the  plasma  and  the  amplitude  of  the  initial  momentvun  in- 
duced. We  would  expect  the  frequency  of  the  plasma  to  vary 
as  the  square  root  of  the  particle  density.  From  the  fol- 
lowing figure  in  which  A'  was  varied,  this  variance  is 
established. 
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PLASMA  frequency (1x10 


One-  and  Tvo-Dlmensional  Runs  with  the  Beam 


An  Initial  time  was  calculated  in  order  to  place  the 
position  of  the  beam  just  in  front  of  the  left  edge  of 
the  material,  and  the  program  was  started.  There  should 
be  a small  initial  displacement  of  the  electrons  toward 
the  beam,  and  then  as  the  beam  passed,  a rapid  acceler- 
ation in  the  beam  direction.  These  directional  displace- 
ments would  be  detected  first  by  an  increase  in  the  mom- 
entum in  the  negative  direction  (beam  sucking  electrons 
toward  itself)  and  then  small  changes  in  the  density  dis- 
tribution. The  internal  electric  field  should  increase 
in  the  restoring  direction  as  the  momentum  starts  to  change. 
These  three  factors  should  allow  us  to  analyze  the  movement 
and  forces  which  act  upon  the  electrons.  The  magnitude  of 
the  electric  field  at  the  various  positions  sho\ild  also 
correspond  to  the  change  in  momentum  as  a fuxiction  of  posi- 
tion. Since  the  beam  is  initially  at  the  left  ( - <»)  , 
the  external  electric  field  should  be  directed  to  the  right 
or  positive.  The  motion  of  the  electrons  should  be  ini- 
tially toward  the  beam  (-)  as  should  the  momentum.  The 
Internal  electric  field  caused  by  this  change  in  momentum 
should  be  a restoring  field  to  the  left  (-). 

All  of  these  projections  are  seen  in  the  initial  few 
time  steps  of  the  calculations;  however,  after  a short 
period,  depending  upon  the  Interval  of  the  step,  the  momen- 
tum starts  to  oscillate  as  described  in  Fig.  12.  The 
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momentum  is  initially  a function  of  the  distance  away  from  | 

the  beam,  but  an  unstable  oscillation  soon  develops  over  | 

time  periods  in  which  the  external  electric  field  does  not 
appreciably  change.  It  would  appear  that  if  a small  amo\ant 
of  momentiun  is  continually  supplied  to  the  system  of  elec- 
trons, an  unstable  oscillation  would  develop.  This  phe- 
nomenon could  not  be  explained,  and  a method  for  damping 
the  oscillation  was  sought. 

It  was  felt  that  the  energy  induced  by  the  external 
field  was  not  being  transferred  properly  into  heat  or  pres-  ^ 

s^lre  without  a small  collislonal  term;  therefore,  the  pres- 
sure term  (p^)  was  analyzed.  It  was  felt  that  if  the  stress 
tensor  was  added  to  this  noncollislonal  form  of  the  pres- 
sure 

Pressure  = p-a^ 

where  p = as  described  in  earlier  sections 
0^  = the  stress  tensor 
stability  might  be  achieved. 

This  stress  tensor  involves  a viscosity  constant  times 
the  derivative  of  the  momentum  with  respect  to  position 

fl(X 

Thus,  the  second  moment  becomes  (1-D) 


- - k + /I  + FN 


f-'. 

i 


I 

: i 


It  was  hoped  that  a value  for  the  viscosity  term  could 
be  picked  which  would  stabilize  these  oscillations  and  allow 
the  hopefully  larger  effects  to  be  seen  as  the  beam  came 
closer  to  the  target. 

Various  ranges  of  V were  tried  without  success.  The 
divergence  of  the  oscillation  could  be  slowed  but  never 
stopped. 

Another  possible  problem  was  encountered  with  the 
solving  routine.  The  ability  of  the  routine  to  establish 
the  proper  8-point  fit  to  the  distribution  and  momentum 
cvirves  was  occasionally  in  doubt.  For  example,  in  one  run 
to  establish  plasma  oscillations,  the  program  diverged  in 
a manner  similar  to  previous  difficulties  with  the  stabil- 
ity in  the  momentum.  A second  run  was  made  over  the  same 
time  region  at  a finer  time  step,  and  the  program  was  stable 
as  predicted.  Thus,  one  would  say  it  co\ild  be  possible  to 
stabilize  the  momentum  in  a run  with  the  beam  present  by 
shortening  the  initial  time  step,  even  though  the  solving 

routine  chose  its  own  working  time  step.  For  several  runs, 

- 8 

the  time  step  was  reduced  to  10  of  a scaled  xmit.  The 
program  wo\ild  run,  but  this  narrow  time  step  would  not  be 
feasible  for  any  computer  routine. 

At  this  point,  it  was  apparent  chat  this  simple  mom- 
ent approach  to  the  analysis  of  the  beam-target  inter- 
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action  needed  to  be  expanded  to  include  either  collisonal 
terms  or  higher  moments  which  would  bring  in  stability 
terms  to  help  dampen  the  divergence.  Artificial  viscosity 
has  been  used  by  several  hydro- code  programs  and  may  also 
be  a candidate  for  this  problem;  however,  a shortage  of 
time  has  prevented  further  analysis. 

Things  to  be  Done 

The  attempt  to  describe  a complex  mechanism  by  a simple 
analytic  tool  did  not  work;  thus,  the  following  two  ap- 
proaches could  be  followed;  One,  patch  up  the  existing 
program  to  include  terms  which  will  eliminate  the  insta- 
bility; or  two,  start  from  the  beginning  with  a more  detailed 
analysis,  including  higher  order  moments  or  even  solutions 
to  the  Vlasov  equations  for  several  particles.  It  is  felt 
that  a patch  of  this  method  would  not  be  advisable  at  this 
point.  The  linear  approach  as  was  done  in  the  collective 
example  should  be  expanded  to  its  maximum  and  also  possibly 
include  a quasi-linear  approach  in  which  the  second  order 
perturbations  are  kept  and  the  higher  orders  discarded. 

This  approach  should  provide  insight  into  not  only  the 
moments  but  to  the  physics  of  the  plasma.  The  plasma 
frequency,  dispersion  relations  and  damping  terms  should 
be  evaluated  as  to  their  importance  before  starting  a pro- 
gram design. 

The  other  approach  would  be  to  examine  solutions  of 
the  Boltzmann-Vlasov  equations  for  particles.  This  approach 


would  probably  not  be  appropriate  to  a master's  thesis 
because  of  the  time  and  details  Involved  in  its  solution. 
However,  it  might  be  the  only  way  to  accurately  model  the 
problem. 
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Appendix  A 


Derivation  of  the  Boltzmann-Vlasov  Equation  From 
Liouvllle's  Equation 


The  equation  of  motion  of  a charged  particle  which  is 
acted  upon  by  an  electromagnetic  field  is  given  by  Newton's 
equation  of  motion  (Ref  15) 

V-£). 

Eq  (1)  is  a means  of  calculating  the  unique  trajectory  of 
the  particles  in  terms  of  their  initial  or  present  position 
and  velocity.  Liouvllle's  equation 


O > 


where  ^ ^ /7^ ^ ^ 

is  a similar  formulation  for  the  dynamics  of  a single 
particle  in  terms  of  the  Hamiltonian  operator,  H.  The 
advantage  of  Hamiltonian  mechanics  over  the  statistical 

nature  of  Newtonian  mechanics  is  that  it  can  more  readily 
be  adapted  to  a many-partlcle  system.  Except  for  very 
trivial  systems  (few  particles),  Newton's  laws  of  motion 
cannot  be  explicitly  solved,  and  thus  in  order  to  proceed 
further,  it  becomes  necessary  to  make  some  mathmatical  ap- 
proximations that  will  produce  a more  manageable  equation 
out  of  Liouvllle's  equation. 

If  the  range  of  forces  within  a plasma  is  fairly  long, 
an  appreciable  force  will  develop  on  a single  particle  if 
all  of  the  other  particles  in  the  plasma. 
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move  in  some  collective  way.  These  so-called  collective 
forces  can  be  expected  to  be  much  more  smoothly  behaved 
f\inctions  of  both  space  and  time  than  the  forces  produced 
by  rare  mutliparticle  collisions.  When  this  theory  is 
carried  to  extremes,  the  argument  leads  us  to  suspect  that 
in  the  lower-order  we  could  replace  all  of  the  Inter-par- 
ticle  forces  by  certain  macroscoplcally-smoothed  electro- 
magnetic fields.  According  to  this  point  of  view,  every 
particle  moves  completely  Independently  of  all  other  par- 
ticles but  follows  its  predictable  trajectory  in  these 
smoother  fields.  Physically,  it  is  evident  that  the  choice 
of  the  fields  must  consist  of  those  fields  that  are  pro- 
duced by  all  of  the  particles  in  the  plasma,  and  these 
fields  are  referred  to  as  self-consistent.  These  fields 
are  obtained  by  solving  Maxwell's  equations.  In  this 
section,  we  will  develop  the  collisionless  Boltzmann 
equation  (Boltzmann- Vlasov)  from  the  many  particle  Liou- 
ville  equation  with  the  aid  of  this  Independent  trajectory 
assumption  and  the  fields  developed  from  Maxwell's  equa- 
tions . 

For  reasons  of  simplicity,  let  us  consider  a gas  of 
identical  particles  and  assume  that  each  particle  follows 
an  independent  trajectory,  or  equivalently,  the  N particle 
distribution  function,  f (t)  , may  be  expressed  in  the  form 

- ir  ^ Cn,  Pi  rt) 

t-  \ 
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where  f Is  the  single-particle  distribution  function. 
Liouville’s  equation  for  the  dynamics  of  a many-particle 


system  is 

IfU 

dt 


*/  j ^ fj 

+ 1 - 1 

J*' 


to,  = O 


(a-0 


where  H = the  Hamiltonian  operator  and  the  J refers  to  the 
Jth  particle 


Hr,  -- 


olH 

drj 


^ex  ■ ■= 


UH 

<LPx  * 


For  most  physical  situations  of  interest,  one  expects 
only  a very  small  number  of  inter-particle  correlations  to 
play  an  Important  role.  Thus,  it  is  convenient  and  pract- 
ically mandatory  that  fti)  be  reduced  to  distribution 
fvinctions  involving  only  small  numbers  of  particles  by 
averaging  out  the  coordinates  and  momenta  of  the  remaining 
large  collection  of  particles. 

Since  we  have  restricted  to  only  one  species, 

V 

r c<)  will  be  symmetrical  with  respect  to  exchange  of  the 
coordinates  and  momenta  of  any  two  particles.  The  s-part- 
Icle  distribution  fvinction  F\4.)  (5<<n)  Is  obtained  by 

averaging  over  the  coordinates  and  momenta  of  the  remaining 
/y/-s  particles,  thus 

TT  ( ^^3 


s 

where  TT  ( dp  - ) 


is  the  joint 
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probability  at  time  (t)  of  finding  a particular  collection 

of  s particles  located  at  ri,  ra,  rs  ....rg  in  the  range 

drx,  dra  ....dr  respectively,  and  with  momentum  p,,  Pa  .. 

p in  the  respective  ranges  dp,,  dpa  ....dp  regardless 
s s 

of  the  states  of  motion  of  the  remaining  n-s  particles. 

The  Hamiltonian  operator  is  given  by  the  equation 
(Ref  15) 


^ ^ (r.  >-t)  -h 


A/ 


i-*7 


■f- 


(fi- 


where  A and  <J)  are  the  potentials  associated  with  the  ex- 
ternal E and  B fields . 

It  will  become  more  convenient  to  write  the  Hamil- 
tonian in  terms  of  two  components.  The  first  component 
h^  depends  only  upon  the  coordinates  and  momenta  of  the 
single-particle  distribution,  and  the  second  com- 

ponent, h^,  depends  upon  the  two-particle  distribution, 

/ i \ 

F'  % Extracting  these  components  from  Eq  (3)  yields 


-f- 


7?  • rTj  ^ 
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we  can  write 


From  Eqs  (4)  and  (5)> 

H-t  K + '4  h^j 


then 


M - AH  - ^ ^ 4-  ' y AV.k  I ’C  </ 


but,  since  the  particles  are  indistinguishable 


AR  - -h  f 

^ dl^  -t^a: 


r,A 


Also,  note  from  our  definition  of  Fc^;  and  the  chain  rule 
that 


A£m)  - Mx.  77  i.-  . 
■*G  ,Vr 


f/l-7) 


Combining  Eqs  (6)  and  (7)  yields 


AJ±  r -/- 

oCCj 


«4x  T^/(  \j, 


(A-e) 


By  similar  techniques,  it  can  be  shown  that 


i.  H J,  Fto  . / «( <'’ 


+ y-  Jdx  jT 

I r.  A I A ^ 


Substituting  Eqs  (8), (9)  and 


f^Ttf 
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into  Liouville's  equation  results  in 


“ ^ ila  |,f'  ‘0^  t '■ 


Let  us  now  look  at  each  term  of  Eq  (10)  ctnd  integrate  out 
the  one  body  term  f * . Later  we  will  examine  the  two  body 
terms . 

Integrating  the  first  term  of  Eq  (10) 


I f-i  ^ 7T  }•,  + ^ ]f  f ,• ) 

''  I J-ri.  /Vr  dlt  t^l  /♦ 

'(A-iO 

The  integral  reduces  if  we  bring  out  the  summation.  For  each 
term  then  we  have 


a A 


or  J.P,  f-T 


but 


j'  ^tx  f = 1 


for  each  term  of  the  sum  and 


(1)  = o 


Thus,  the  only  remaining  term  from  Eq  (11)  is 


j f.-  ^ 


t'¥T 


(A-tz) 


since  each  of  the  product  terms  equals  1,  the  total  product 


= 1. 


The  second  term  of  Eq  (10) 


after  integrating  out  the  J^l  terms  yields 


plus 


I f dr,  dri  if  fi''  ( TJ'^hj  f <• 

s*t-  J K-*T  <kTj  \ /^? iJ  J 

Since  the  integral  over  each  f<-  = 1 the  summation  is  written 

(N-i)\c^rU:P  ^ or  df^. 

Integrating  by  parts  leaves 

- fv-,)  / </a  «>.  /. 


where  J = a is  an  arbitrary  substitution  for  any  one  of  a 
set  of  (/V-1)  values. 

Doing  the  same  for  the  third  term  of  Eq  (10) 


l_  /A}L  iR]fj  f.- 

ir,  [ ir,  ] 


yields 


4-  f dr^  d^*.h  • (A-<r) 

J Jin  Jin.  ' 


I 


Combining  Eqs  (13) > (1^)  and  (15)  yields 

— Ah.  AL  4 Tfce  TUJO  %•*«!  -ic» 

Jii  An  An  An  An 


(A-n^) 


Returning  to  Eq  (10)  and  looking  at  the  two  body  terms 

k WG  ' LhK  [ An  <hrTj}}^ 

-TTT  ^ 

rewriting  the  TT  f / as  f k + TT  f < 
and  factoring  out  the  J = 1 term 


Where  7 

L AO  d -pT  Q-  7 

a Poisson  bracket. 


combine  to  form 


Integrating  Eq  (17)  over  all  space  and  momenta  except  1,- 

£ f J ^ 

I j 

all  of  the  other  terms , J <<^/<}  /<  * -2  i r . 

Looking  at  the  first  term  of  Eq  (18) 

f f [ct/,  c/n,  Ah.  - ihl  u 

Jr.hjJ  \j.r,  -TtT- ir,J  ' 


Integration  by  parts  with  respect  to  dp  in  the  first  term 
and  dr  in  the  second  yields,  i.e,  the  first  term 
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then  the  integral  equals 

_o 


Let 


c/«- 


< r _ £J^1  r 

^ J.oo  J 


V-  ' 4^‘/» 

J.r 


(A -1 9) 


the  second  term  yields 


^ Ji^  l^j-f. 
d-p  dr 


so  the  first  term  of  Eq  (l3)  disappears  and  Eq  (12)  reduces 
to 

jd 

( t ^ j t ^ iT*^ 

(a-  19) 


i j f_ 


Since  h*  is  a function  of  (p*,  p^,  r*,  r^)  all  terms  of  the 
sum  will  be  the  same.  Therefore,  Eq  (19)  reduces  to 

(N-i)  { T-dilli  ^ >(  f'l  1 

‘ J * LT^  'Z7T- rrrJ 


where  K - 2 is  a dummy  selection  or  rewriting 

Tr.  '■>  / Qik  -^L  (^">1  *•'  fj z^r 

where  (a/-i)  A»t  is  the  potential 

due  to  all  other  (n-1)  charges  on  particle  1,  Therefore, 
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r 


z (A/~i)f<iK  c/Cj,  4*.  (a-io) 


Combining  all  the  terms  , Liouville's  equation  becomes 


jJj  - A}i  AL  - ill  iL  + ^ ^ 

d ‘i-  d."/)  dfi  d fi  d(~i  d'Fi  ^ JIO 


If  we  let 


h'”  = A'  * V 


< ,/l 


then 


^ c(r  e/f  h'- (rr'/ f>^') 

= *i'  (r,^;  + v(r)  . 

Excluding  magnetic  effects, 


where  (|)  is  an  external  potential  so  that 


-sU!‘^  - -£. 


«A  h 


J.r 


11  4.  dl. 
dr  dr 


and 


AL 

di 


- f A'"  X ? - a - 

■ C ^ ^ ■ L I^r  d-P  di*  dr) 
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but 


then 


_a1  . fjA  ^ 4j^]aL  - Al 

di  'Wr  Jr ) ii>  AA  dr 


^ r - V • F”  4 

AA  dr  V j'  T 


[-^  + -jf]  - (-  /^»rc€^fx^)-  /rcrc<.(XAji))^ 


therefore , 


-f  V • 7^  -k  J •>•'<.<.  • ^ j s o , 

A 


The  internal  force  A-**  = ■"  V { 4 C'‘.P).»  ^ , 

■J  I 


letting 


r - V r J^r  /^(r)c^  = e £T  . 

J v/rr,  /r-r*/ 

Therefore,  Liouville's  equation  becomes 

^ + v*7f  +ec^*rp4=®- 

If  magnetic  effects  are  included,  h becomes 


. I 


(p--g./?y  + 
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i i 

i! 


^ where  is  the  magnetic  potential  and  Q - ('YX/7) 

it 

Through  a similar  process  of  taking  the  derivation  of 
t « with  respect  to  p and  r and  substituting  into  Eq  (21) 

Llouville's  equation  becomes 

^ + V • 7^  + e [e  + (VX0)J  • dl  z o 


the  Boltzmann- Vlasov  equation,  where  the  symbol  V takes  oh 
the  form  V=  (P-e/7)  instead  of  V = -H.  . 
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The  External  Electric  Field 
The  external  electric  field,  due  to  the  presence  of  a 
IK  amp  beam  of  10  GeV  protons  chopped  into  pulses  of 
10“ ® seconds,  has  been  modeled  as  a point  charge  of  10“ ® 
coulombs  moving  at  a velocity  of  2.98  x 10®”^/s.  The  assumed 
radius  of  interaction  has  been  selected  as  .001  m or  in 
other  words,  the  radius  of  the  beam  is  assumed  to  be  on  the 
order  of  .001  m,  but  since  we  are  modeling  it  as  a point 
charge,  the  area  of  greatest  concern  will  be  the  cylindrical 
area  .001  m from  the  axis  of  motion  of  the  point  charge. 

In  this  section,  we  will  develop  the  equations  for  the  axial 
and  radial  components  of  the  relatlvisitlc  point  charge. 
These  equations  will  be  used  in  the  calculation  of  the 
Lorentzian  force  felt  by  the  electrons  in  the  plasma. 

The  electric  potential,  (|>,  due  to  a point  charge  at  a 
distance  (r^i,  is  given  by  

where  is  a retarded  charge  density  of  a point  charge 
given  by 

z Q ^ (Z-Vi)  ^ (r) 


where  the  point  charge  is  represented  by  the  delta  function 
in  radial  and  axial  dimensions,  V = the  velocity  of  the  beam 
and  t = time. 


1 
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V. 


c =s  the  speed  of  light 


I r-ri|  = the  distance  from  the  point  of 

interest  to  the  beam  as  shown  in 

Fig.  19. 


0(R,t)  = Xcc)  Xrz  - VC-t  - 


and  Ir^l  = + (x-z)^  • 

Let  us  define 

rC2)  = 2 - V-d  - /?  fx-z}*- 


then 

^ 

for 

/r(t)  r O = 

= <? 


thus, 


rzo  - / ^ /^Cz-X) 

U Z +(k-z?~  (8-3) 


- £2.  -I-  V'li"  = ^"r''4-  - ZP'XZI 


Rearranging  with  powers  of  Zo  yields 


- £.(2vt  -z/»x)  /V*;  = o . 


If  we  now  let 


= (l-r)  , , C =rv/‘r-/^V"-/9V‘; 


and  use  the  quadratic  equation  to  find  the  roots 
2.  z ^ rvt 


/ ^2 - /3*xj " /»*;  /s 

2-r;-'n  ■ 


- /3^x  ± 

('/-/3V 

In  order  to  determine  which  sign  should  go  with  the  square 
root  term,  let  us  add 

(-Vi  + + Cx-v^j*- ) 

to  both  sides 

^ _V-d  4-/3  V - O = V-f-  ± ^f~(\-P') 

^ (8-^) 

In  the  line  as  (vt-X)  goes  to  zero,  the  right-hand  side 
goes  to 

O z ± + R /(/-rV  •*-  /®*’(/-/®V 

therefore,  it  is  clear  that  the  sign  must  be  negative.  Then 


- fi'x-  p *(v-t-x)" 


I - (9 


(a-i) 


From  Eq  (3) 


o(  Ri») 

I z 


V +('X-zJ)^  + /^(Zo-.x) 

^.(X-Z.)*-  ((3-8) 


But,  Eq  (5)  says 


R'"  4-(X-Zo)^  = 


thus . 


F g.)  . *■ 

^ Zo  ~^± 

fl 


Zo-M-i  ^ P'U.-K 
(Z.  -V^J 


(b-*») 


Returning  to  Eq  (7)  and  subtracting  (vt)  from  both  sides 

I - 

and  subtracting  (x)  from  both  sides  of  Eq  (7) 


(io  - 


I - 


Substituting  Eqs  (10)  and  (11)  into  Eq  (9)  yields 


/ F (g.)  - 


Z^(Vi-x)  - (Vt-x)" 
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Combining  Eqs  (2)  and  (3)  and  noting  that 


-(Zo  -vt)  = I r'-Tl  1 
6 


CS  (R.o  = ^ f t 

^ v-£)  / , ^ p> 

\ ^o 


7 


j 

I 


(8-1^) 


As  expected,  Ez  reduces  to  the  familiar  form 


*"  " VTTf. 

as  6 -»•  0,  and  the  area  of  the  interaction  goes  to  a point 
(R^o) , In  the  one-dimensional  model,  the  radius  of  inter- 
action will  he  fixed  at  (.001)  even  though  the  r dimension 
is  nonexistent,  and  in  the  2-D  model  it  will  be  allowed 
to  vary. 
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Appendix  C 


Because  of  the  large  number  of  parameters  and  scaling 
factors  used  in  the  analytic  solution  and  the  equations 
which  are  used  in  the  computer  model,  this  section  will  go 
through  the  rationale  used  in  the  selection  and  derivation 
of  these  values. 

The  following  constants  have  been  used  in  evaluating 
the  one-dimensional  Vlasov  equation; 

1)  Xo  , the  distance  scale  parameter  where  X = Xo 
Assume  target  to  be  1 meter  thick  and  extends  symmetrically 
over  17  cells  with  values  of  u =-•/<»  ^ +%  in  steps  of 
.111,  thus 

2)  , the  number  of  ele ctrons //«?■*  - momentum 

where  N = mimber  of  electrons/^*  in  the  plasma.  Assvime 
the  plasma  electrons  are  the  3 outer  shell  electrons  of  the 


aluminum  atom  then 


X /.8  X 


It  should  be  noted  that  if  energy  greater  than  a few  hundred 
is  deposited  in  the  plasma  per  aluminum  atom  more  than  the 
three  outer  electrons  will  be  in  the  Conductance  band.  For 
the  first  estimate,  only  the  outer  three  were  chosen  but 
additional  analysis  of  future  data  may  require  this  value  to 
change . 

Po  is  a free  parameter  and  will  be  discussed  later. 
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■iggiPWiNwi^PP 


3)  Eo  > the  electric  field  scale  factor  which  will  be 
discussed  in  Appendix  D and  serves  to  condense  the  computer 
operations well  note 


Pq  — 1 . e.  E o _ f 

M y.  ' p«  ■ 


and 


To  e.  P/  ^ = j_ 

^e. 


Solving  these  equations  for  ro  , p„  and  £©  yields 


T;  r 7.ZS- 


P.  I 3 V 0 X /O  - ^ 

Eo  T 3.0  X /«*"  V//m 

4)  C?  , the  charge  of  the  pulse.  The  beam  consists  of 


a IK  amp  proton  beam  chopped  in  10"’  sec.  pulses.  Thus  each 

V-  6 


pulse  has  a charge  of  <3  = 1 x 10  coul. 

5)  fi  ) the  relativistic  constant.  The  protons  are 
rated  at  10  GeV  or  10  “ev.  From  relativity  theory  we  know 


fAo 


where  T = kinetic  energy  of  particle,  M»=  rest  mass  of  pro- 
ton. Therefore  A = -995  and  the  velocity  of  the  beam  equals 

V-  efi  - 2.  70  i 
The  two-dimensional  solution 


1)  l?t> , the  radial  distance  scale  parameter.  If  we 

.-3 


assume  that  the  diameter  of  the  actual  beam  in  10  ' meters 
instead  of  a point  charge,  the  area  of  the  target  which 
will  be  of  greatest  interest  will  be  within  a radius  of 
5 X 10"’  meters.  If  we  have  nine  points  of  observation 
in  the  radial  direction,  then  8 of  those  should  be  in  the 


5 X 10" ’m  radius.  Therefore: 


r\  U-T 

OiO005  = 7_7Z>' 


Let  uy=  V7  then  #?o  = 6.25  x 10 


- s 


1 
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2)  Eo  , T. 


2)  Eo  > T"-  » Po  > Bo  , the  electric  field,  time,  mo- 
mentum and  magnetic  field  parameters  respectively.  As  will 
be  developed  in  the  two-dimensional  method  section  (Appen- 
dix F),  the  nine  constants  of  the  solving  equations  are: 


Ho 


. PoT, 


—7^ — 

/‘1  - e Po  fJo  T* 


C £■  - !* 

c.  R. 

C>.  f?o 


By  making  the  following  assumptions: 

^ t - Cjr-  C7  = Ce  = 1 

The  nine  coefficients  reduce  to 

G = / <r*  r Po  /z-o 


res  / 


C7  r i 


to  To 


'*■  8,  To 
£ • to 

Eo  r, 
&o  2; 


The  values  of  the  scaling  parameters  thus  lead  to: 

- 2.0  7,0" 

Qo  ~ ^o/c  ~ ^‘7Q 

Pc  - Aa  - 2 73  X «'<>■**■  **V* 

a* 

To  3 A\  Po  _ 2. « 8 x/o-'^  s . 


3)  All  other  parameters  remained  the  same  as  the  one- 
dimensional  model. 
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Appendix  D 


r 


r 

L. 


Moments  for  the  One- dimensional  Model 
The  taking  of  moments  involves  the  integration  out  of 
the  momentum  compoents  (phase-space  components)  in  the 
particle  distribution  function  and  is  given  by 

J J a>) 

where  ai  = corresponds  to  the  1st,  2nd,  etc. 

moments . 

The  first  three  moments  of  the  Vlasov  equation 


o/fcx.P)  o-  P •(]■(<••*)  X p)_ 

—Tt ^ -Tx  ^ ^ -3rr-  - 


(p'f) 


in  one-dimension  will  be  taken.  As  stated  earlier,  the 
force  F,  is  the  one -dimensional  Lorentzlan  force 

to- 2.) 

where  £"  is  the  total  electric  field  due  to  external  charges 
and  currents  and  the  internal  field  due  to  displacement  of 
any  charges  in  the  plasma. 

Once  the  moments  have  been  taken,  the  electric  field 
components  of  the  force  will  be  examined,  and  the  scaling 
parameters  discussed  earlier,  will  be  introduced  into  the 
equation. 


The  First  Moment,  m = o 

Multiplying  through  Eq  (1)  by  dp  and  integrating  over 
all  momentum  ( ± oo  ) yields 

^ f-f(x,p)c/p  = -7k  [si(f  ^ ^ ‘^£[- 


IQlf 


where  the  derivatives  with  respect  to  time  and  space  are 
independent  of  momentum  andhave  been  taken  out  of  the  in- 
tegral. Let  us  now  define 

N - f f fx.  P)  p 

where  A/  is  the  particle  density  and  is  obtained  by  inte- 
grating the  distribution  function,  jcji.p')  , over  all  momentum 

A/<  P>  r J p j:  (X,P) 

where  p is  the  expected  value  of  p and  is  defined  as  the 
average  value,  p , plus  any  small  perturbations,  "p  . 

P + P 

So  TA/ir 

<p>  - <p^li>  = <p>  + <ff> 
but  by  definition  ^p>=  p tfcoj  <'p>  = o . 

Eq  (3)  thus  becomes 

(-> 

The  second  term  on  the  right-hand  side  is  a perfect 
differential,  and  is  set  to  zero,  if  we  assoime  that  there 
are  no  particles  at  too  . Eq  (4)  thus  reduces  to 

All  : ^ ± r £l  + p aULl 

Jii  ~ ^ JLX  ^ ^ X J ' 

which  agrees  with  Seshadri's  equation  Eq  8.9  (Ref  l4)  for 
a one -dimensional  model.  Dividing  through  by  /V  , we  obtain 
the  first  moment 

<0-^^ 
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The  Second  Moment,  m = 1 

The  second  moment  is  obtained  by  multiplying  Eq  (1) 
by  pdp  and  integrating  over  all  momentum  to  yields 

-^1  Pf^p  . -.,^1  P‘j  • (0--T) 


Eq  (7)  reduces  to 


(0-8) 


I 

By  a similar  analysis  as  in  the  first  moment,  pf  I ^ - O 
< P*>  r < p'  + 2 P P -f-  'p*')  - P *■  • 


Thus,  Eq  (8)  reduces  to 


^ p^-h  <T))J  + p/V.  (o.<^) 


Let  us  now  define  as  the  kinetic  energy  term  of  the 
plasma.  Since  P =^"^the  term  may  be  written  as 
2^*"  = //(mv‘)  which  has  units  of  kinetic  energy  per  unit 
volume  and  a/<  is  a kinetic  pressure  term. 


If  i.  P>  is  written  in  vector  notation  as  < p > = 

^ Pf  Pt}  > can  be  shown  that  only  the  following  terms 
are  nonzero  in  three  dimensions 

< P > ^ < Px  Px>  f <PyP,'>  ^ < PiPi) 

j/j  < _ K-r 


fA^ 


in  1- Dimens ion 
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where  K is  the  Boltzmann  constant  and  T is  the  uniform 
temperature.  The  static  pressure  is  given  by  KT  thus 


M. 


<f> 


/v\ 


P. 


Thus  Eq  (9)  may  be  written 

= -jVC'VP 


r) 


(D'/o) 


which  agrees  with  Eq  (II-3)  of  La  4700.  (Ref  17)*  For  our 
purpose,  however,  let  us  keep  the  notation  of  'p*'  and  assume 
it  to  be  a modified  pressure.  Thus,  Eq  (10)  may  be  written 

expanding  Eq  (11)  and  substituting  into  the  expansion  of 
Eq  (5) 


(D-IX) 


- + A'4t  + 


Canceling  out  the  A/  yields 


j/JB  - pr  ~ J_  r p ^ 4-  e/  P I n ol  P 

4-t  - ^ M ^ XT  ^ P Xx-J 


(d-  »3) 


the  equation  for  the  second  moment. 


smmaassm 


The  Third  Moment 

Similarly  as  in  the  first  two  moments,  Eq  (1)  is  multi- 
plied through  by  p^dp  and  integrated  to  yield 

= - wljhti  - //  ‘^P- 

After  integration  we  obtain 
^(A/Cr+P*)J  r -X  ^£y,/<p3;J  + £ ^ 

But  <p^>  = <^*+3pp  + p*>  when  the  <p>  terms  are  set  to  zero. 
Eq  (Ilf)  then  becomes 

( A/C<P^i-  P‘))  ~ --Jk  -hi  + ^ -hZ^'^P  ■ 

(D-ir) 

In  order  to  close  the  set  of  equations,  let  us  set  all 
p’  terms  equal  to  zero.  In  actuality  this  term.  A/  , is 

/v\ 

an  energy-flux  term.  Eq  (1?)  then  becomes 

-f^[N(p\<p>)]  - ~ ^ -t  z^NP  . 

(o-/fc) 

Substituting  Eq  (11)  and  Eq  (5)  into  Eq  (16)  and  collecting 
terms  yields 

if--  (o-n) 
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The  third  moment,  which  agrees  with  Seshadri's  equation 
(10.31)  if  written  in  Eulerian  derivatives  where 


To  quickly  siimmarize,  the  first  three  moments 
m = o 


Let  us  now  look  at  the  force  term,  F.  In  one  dimension, 
term  is  given  by 

F = 'l^Eint  + Eext)- 

The  external  electric  field  is  developed  in  Appendix  fi 
and  will  receive  no  further  discussion  at  this  point,  ex- 
cept to  note  that  it  is  Independent  of  the  plasma. 

The  internal  electric  field  can  be  evaluated  through 
Maxwell's  equations 

(V  XB)  ^ I + €0  ^ (o-ie) 
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where  7 is  the  current  density  and  is  given  by 


I 


J ^ ^ I P ^ P)  • (O  - (<») 


Combining  Eqs  (18)  and  (19)  and  noting  that  (V  x B)  = 0 
in  one  dimension 


■ir  - 

After  integrating  over  all  momentum  Eq  (20)  becomes 


AL 

It 


a 

£»A^ 


/\/p 


(^1,-21) 


where  q = -e  the  charge  of  an  electron. 

The  Lorentzian  force  may  then  be  written  as  a sxun  of  the 
two  electric  fields  after  a differential  equation  solving 
routine  has  been  used  to  calculate  Eint» 

At  this  point,  we  have  four  linear  differential 
equations,  Eqs  (6),  (13) » (17)  and  (21),  which  can  be 
solved  by  a standard  solving  routine  and  thus  stepped  in 
time. 

As  was  discussed  in  a previous  section,  it  is  necessary 
to  scale  the  parameters  of  the  four  equations.  The  seal- 


Pa 


dP-po  foSK  (V) 

E = (’f;l 


f - 

E =£.£.* 

With  these  substitutions,  our  four  equations  change  to 


A/  _ 


M 


C^5k(S)f Caik(v)M  4-  ^Jn^l 

L-  •(  J 


S/VK  CV) 


^v'  - To  cE„  c<3SAfs)E*-  TrJ?  cdsk^s) 


A^Xo 


£ S/A^k  (✓‘J  + CoSKCvV^  + S'»^(.'^)  Coik(v) 


(b-2^3) 


S.a/KCv'J^'  r ^ 1^2  C«sk(v/)  ^ + s///krv)cosK(/V  ^ J 


(D-2«J) 


4 c - t; e Pg  Vo  r (2V)  c/v 

J s " J 


(D-^^'> 


I 


If  we  now  let 


To  Po  ~ 1 

M 


To  e,  - 1 

Po 


T,  a P^A/c  _ /I 

z t,  A1  E. 


these  three  equations  can  be  solved  in  terms  of  one  para- 
meter, Eo 

Eo  = Xo  e No 
2eo 

Call  Yi  = In  N 

Ya  sinh(v) 

Y3  = slnh  (v^) 

Y4  = Eint  • 

Then, 

Eqs  (22),  (23),  (24)  and  (25)  reduce  to 

r - R [ y,  + ^ + y.  4 4 E„.) 

(D-il) 

a = ^ + >4 


(d  - ^8) 


Appendix  E 

Derivation  of  Moment  Equations 
In  Cylindrical  Coordinates 
In  this  section,  a general  discussion  of  the  metric 
tensor  and  other  properties  of  tensors  necessary  to  find 
the  moments  of  the  Vlasov  equation  in  2-D  cylindrical  co- 
ordinates will  be  discussed.  Once  the  basic  groundwork 
has  been  established,  the  first  three  moments  will  be  taken. 
Maxwell's  equations  will  also  be  integrated  into  the  equa- 
tions as  required  by  the  Lorentzian  force  term.  These 
equations  will  then  be  simplified  to  meet  the  assumptions 
of  our  model. 

The  Metric  Tensor 


The  distance  between  points  P and  Q is  ds . To  evalu- 
ate this  distance,  we  take  the  square  root  of  the  sum  of  the 


i 


11 


If 


H: 


M- 


I 
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r 


squares  of  the  three  increments  dr‘,  dx*  and  dx* , 
(ds)*  = (dxM^  + (dxM^  + (dxM'. 

If  X has  the  form 


then 


= f^(z‘,2Sz‘)  1 = 1,2,3 


41 

<k  * 


where  summation  notation  has  been  used  and  s represents  a 
s um  s~  1^2j«**etc«  but 

(ds)^  = dx^dx^ 
so 

(e-i') 

where  the  dummy  index  s has  been  changed  to  avoid  confu- 
sion. Eq  (1)  represents  a homogeneous  quadratic  form  in 
the  variable  dz^.  Let 

b,.  - (Jl^.  at  ) 

~ \ dZ‘  / T = I.  2.  3 

then  (ds)^  has  the  form 
(ds)*  = bj^j  dz^  dz*^ 

where  bij  for  a three-dimensional  space  is 

(biv  bi2  bis' 
bai  baa  has 
bs  1 bs2  bs  s< 

The  expanded  form  of  the  quadratic  equation,  Eq  (2),  is 
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(c/s)  - h.io/to/i'  + biz.  c/g'c/z*- 

4-  bzf  d.1^  ei.%'  4 • • • ' • 


b.3  o/z  o(z^ 


+ bfi  xs/e* 


Rearrangement  yields 

U/s)  - bn(c(z')^  4-bzzCa/l')  f bs3  (c/z*)  t(bil-i-hzi) 

yo/z'e/i*  f (hi3  + bSi)  a/z'Wz'^  -f(bz3  -^b3z)  c/^o/z^ 

or  in  a general  form 

(e(  %Y  - e/z'  . (e-^ 

Let  us  define 


~ (bir  -t  bjt  ) (E  - 

Z 


hence , 

(Vs)^  ^ V'a' 

(Note;  g,’j.  = g-j.). 

ZIj  is  a covariant  second-order  tensor  and  is  called  the 
metric  tensor.  A covariant  tensor  transforms  a tensor  in 
one  system  to  a tensor  in  a second  system  via 


The  Metric  Tensor  In  Cylindrical  Coordinates 

The  Cartesion  coordinates  , x-*,  x*  are  related  to 
the  cylindrical  coordinates  r,  (j),  z by 

x'  - r CC75  © r i?'  cos 

X*-  r e - 


Employing  the  general  equation  (Eq  (^-))  yields 


O/,  - A(’Tco'se>'i  ‘l(^C£>se>')  J.(  UCrsliJe) 

^ ^ ' Ji  z*  ^ ^ "z ''  ■ — 12^ 


If  ( = J=1 


4 i 


r '^Zr  - c/r' 


then 


&'  = / =■  J/ 


i ' T ^ z. 


If  / - J-r  3 


= o o + 1 ^ 1 ' ^3: 


all  other  possible  combinations  of  i and  J result  in  g^.j. 

Th«.  3u  = (t  r‘=j 
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(c-s) 


This  is  the  metric  tensor  in  cylindrical  coordinates. 
Thus,  the  distance  between  two  points  in  cylindrical  coor- 
dinates is  given  by  (ds)^  = ^ or 

/,  O O V 

(c(s)  - fo  d^i  . 

\o  o I / 


Chris toff el  Symbols 

The  partial  derivative  of  a vector  A with  respect  to 
is  given  by 


^ ^ s.  + jA  /,'• 


(e-  O 


where  A = A a,  and  ai  is  a unit  vector  of  the  space  to  be  used. 
Let  02s  define 


Where  ^/r]  is  called  a Chris toffel  symbol  of  the  second 
kind.  Substituting  Eq  (7)  back  into  Eq  (6)  yields 


UK' 


^ • 


In  other  words,  the  Christoff el  symbol  serves  to  take 
the  derivative  of  the  space  in  which  the  vector  is  defined. 
Thus,  the  covariant  derivative  of  a covariant  tensor  is 
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defined  as 


R ',r  Z 


(C-  8) 


where  j means  the  covariant  derivative  with  respect  to 
J , and  a covariant  tensor  is  defined  by 


8. 


where  is  defined  in  the  y- coordinate  system  and  A in 
the  X-coordinate  system  and  the  coviariant  derivative  of  a 
contravarlant  vector  is  given  by 


where  the  superscript  means  contravarlant,  and  subscript 
means  covariant,  and  the  contravarlant  vector  is  defined  by 


a'  = /7‘ 


If  the  tensor  is  of  rank  two  (i,j)  Instead  of  rank  one 
(a  vector),  Eq  (9),  expands  to 


In  Ref  18,  the  Chris toff el  symbols  for  cylindrical 
coordinates  are  derived  and  will  be  only  stated  here  as 

[ = { 1^2) 

{ 2^2  } = ■ 


Equipped  with  the  Chris toff el  symbols  and  the  definitions 
of  the  derivatives  of  tensors,  Eqs  (8),  (9)>  and  (10), 
we  are  ready  to  take  the  moments  of  the  Vlasov  equation. 


The  Moments  to  the  Boltzmann-Vlasov  Equation 
The  Vlasov  equation  is  given  by 


(kt 


/e  - 


where  f = the  distribution  function 

- the  moment\im  tensor  in  the  <<  direction 
^ = the  mass  of  the  electron 
F*  - the  force  tensor  in  the  <<  direction. 

If  we  define 

yf  (/)  = / ^(XjP)  (e-'3) 

where  '>t.3  and  is  the  value  of  the 

average  momentxim  as  was  discussed  in  Appendix  D. 

Then  after  multiplying  Eq  (12)  by  and  integrating  over  all 
momentum  we  obtain 

where  r'is  a function  of  X and  t usually  thought  of  as  the 
Lorentzian  force^  f + (Vxfl))and  is  independent  of  the 
integration. 
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The  last  term  on  the  right-hand  side  of  Eq  (1?)  vanishes 
since  it  is  a perfect  differential,  and  the  value  of  ^ at  the 
limits  is  assumed  to  be  zero. 

Substituting  Eqs  (13)  and  (14)  into  Eq  (15) 

/f  (x)  T - -L  /f  ^ 

(c-  /fc) 

which  is  in  agreement  with  Eq  (8.9),  the  equation  of  con- 
tinuity, of  Ref  14.  Expanding  the  covariant  derivatives  via 
Eq  (9)  yields 


(X.  = r,  ((),  z = 1,  2,  3 respectively  or  x,  y,  z if  in 
Cartesian  coordinate,  only  the  values  of  the  Christoff el  sym- 
bols change.  The  assumption  is  now  made  that  the  momentum  in 
the  <J)  direction  will  always  be  nonexistent  due  to  symme- 
trical modeling;  thus 


0 

and  Eq  (16)  becomes 


Jit 


(£-  '^) 


where 


Define 


^ (T -/?")( 


where  (T-/9)  Is  a small  perturbation  term,  = (/f  - joj,  and 


/f  is  the  modified  pressure  term  as  described  in  Appendix 


D. 

Upon  expanding  the  integral,  the  equation  reduces  to 

/i/C"  = f 

Multiplying  Eq  (12)  by and  integrating  yields 

^ f 

or  after  substituting 

- -sf£/fw4w+  /?w  /?*«•> /f*")  2^ 


y/* 


delta  function  = 1 if  o(  = I?  , 0,  otherwise. Eq  (II-3) 
of  Ref  17 > if  the  force  term  is  added,  agrees  with  the  above. 

Expanding  the  differentiation  and  substituting  Eq  (16) 
yields 

/-  -k  I f/r?  /?’ 

yM  y - y ^ y/’'' [M . 

Thus , for  ^ - r o(  = (*”,  <^,Z) 
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Let  us  now  define 


-- 1 ‘^pfr'^pY-  p.pY  - p- fi  n. 
-pp*Y  - v/s/fV''/’’' 


The  third  moment  is  found  by  multiplying  through  Eq  (12) 
6 ^ 

by  P P and  integrating  to  yield 

-f,  [ rYi  --  z-W/  ^?] 


Substituting  in  definitions  result  in 


--  ^:fY  f.7. 


ifr  .!», 


* f,Y  ^ +P^LF.rJl 


4-  f'^LFF"]- 


(E-zz) 


The  left-hand  side  of  Eq  (22)  changes  to 


f>  i p:\  + r:r.4^  /f'  4 /f'/f  i,?! 


Solving  for  <</^  by  substituting  into  Eq  (22)  the  simpli- 
olt 

fied  values  derived  or  defined  above  yields  after  much  sim- 


plification 


^ P - -nr 


Cl 


cCx*< 


Xt' 


4.x 


(e-z3) 


Now  if  r = B = r 


/f"  - - 7^  r "4^  (''5"  * 
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_J_  « 6 -h  +p/^  ~ ^ fz  ~ -4-  LI°I 

j^r  elZ 


(1^/4  ff  f £^7%/?y 


+ - r r />/% 


+ Lfi^r  + /i'‘/f%/fV'‘)C  -rfl'^J 


+ /fY-ir- ) + f‘\-rP^''-  rps'*'\  + /? 


(E-  2*4) 

In  a similar  manner,  Eq  (23)  must  be  solved  for  /- 

= y = H)  r = £. 

The  total  process  results  in  ten  equations,  but  in  the 
model  we  have  chosen,  there  is  symmetry  in  the  <f>  direction 
initially  and  remains  throughout  the  problem;  then  ^ 

equal  zero.  The  P^*** ^ terms,  however,  do 

not  disappear. 

With  the  above  simplification  and  stopping  with  the 
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+ if'J  * 2 f 


#'7  - 

fi^afj^Pi 

if  7 7 

Cf  - 28) 

it 

^ il  J 

(k-  s.'^) 

1 - 

- _ / r <kJ^/i 

4 

+ /?"  ^/?v'  * 44!j 


(E-3o) 
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i4-  '*• 


( 


A 


+ 


The  only  terms  yet  to  be  defined  are  the  force  terms 
^*"and^*.  The  forces  exerted  on  the  electrons  will  come 
from  both  external  and  internal  electric  and  magnetic  fields. 
For  this  model,  we  have  chosen  the  external  magnetic  field 
to  be  zero  because  of  the  large  difference  in  the  magni- 
tude of  electric  and  magnetic  fields  of  a common  point 
charge. 

The  force  on  an  electron  is  given  by 

F - s{i  Vxb)  (e-«) 

Wk.rt  g - , g 

E external  is  evaluated  in  Appendix  B . The  internal  electric 
and  magnetic  fields  can  be  obtained  from  Maxwell's  equations 
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(E  ■ 33") 


- Vx  E 


^ VXB  -■£  • 

.O 

Since,  as  we  have  stated  earlier,  only  the  ^ component  of  B 
is  nonzero  and  only  the  r and  z components  of  E are  non- 
vanishing, then 

- -P  - AL±1  (e-3s) 

Jit  - L <kz  JLC  J 

where  cylindrical  coordinates  is  given  as 

(r-34) 

The  time  derivatives  of  the  electric  fields  are 

^ Er  - _ _L  AlsL 

(kt  to  ^»io  tk  ^ 


c(Q 


# -- 


- J_ 
A 


thus. 


c/  Lr  _ ( P AMsL  f s 


where  C = speed  of  light 
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^ - 'rr  ^ 


: /f'/?  4 4^). 


(E-  38) 


The  differential  equation  solving  routine  used  to  solve 
for  the  moments  can  be  used  to  obtain  solutions  for  £"r  , E 2 
and  Equipped  with  these  solutions,  the  forces  on  the 

electrons  can  be  evaluated.  The  forces  in  the  radial  direc- 
tion are  given  by 

Jr  ^ -e  (Er  -h 

- -e(Er  - 

where  = Er*'^  + tr 

The  forces  in  the  2 direction  are  given  by 

-e(^i  r -e  (E.  4-^'’s^j 


In  s\immary,  we  have  ten  partial  differential  equations 
with  respect  to  time  which  when  given  initial  conditions  can 
be  stepped  in  time  to  yield  solutions  at  later  times.  The 
ten  equations  are  Eqs  (25-31)  and  Eqs  (3?)»  (37)  and  (38). 
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Scaling  the  Equations 

As  In  the  one-dimenslonal  model,  It  becomes  necessary 
to  scale  the  various  parameters  In  order  to  reduce  the 
large  magnitude  of  the  numerical  values.  As  In  the  one- 
dlmenslonal  case, 

T - To  r 

^ p.  V 

X - XoX*  or  r ’ R«  ur  ^ 2 . u-  ct 
B = E,/»7 
e * B.  ^ 

f “ JZJ  f Ktrt  ^ ' 

Thus,  Eq  (12)  may  be  written 

if  - - ('''^'>)  3T- 

(e-W 

In  a similar  manner,  the  equations  for  the  time  de- 
rivatives at  the  electric  and  magnetic  fields  are  given  by 


MEof.  J 


(c  - 


(e.uO 

since  we  may  not  wish  to  work  with  the  same  scale  para- 
meter In  the  r and  z direction,  we  can  expand  the  above 
equations  to  have  the  same  form  except  the  terms  will 
have  two  components,  and  and  corresponding  values  of 
Ro  and  2.  . 
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k , 


1 


Let  us  define 


r - Ml. 

~ AA  Ro 


~T. — 


Cy  = 


gB»/» 


Ci- 


B.  /. 

£..  1?* 


Cfc  - 


c*-Bor» 

to 


^ - c Po  M>  7» 

M£.  £o 


- 

® " BoR 


C«» 


Eo  To 

8o  2.  ' 


These  nine  equations  have  only  four  unknowns  if  we 
pick  Po , Zo  to  be  related  to  the  target  dimensions  and  can 
thus  be  solved  to  yield  the  values  given  in  Appendix  C. 

Eq  (39)  thus  reduces  to 


and  the  electric  and  magnetic  fields  to 

-f 


(t  - 


Eqs  (4-2),  (43)  and  (44)  have  the  same  form  as  their  un- 
sealed equations  used  in  previous  sections  of  this  appen- 
dlX|  except  for  the  constant  which  can  be  carried  through 
without  modifying  previous  work. 

As  was  the  case  for  the  one-dimensional  model,  the 
true  form  of/,  \/,  X , / is  given  by 

P - Po  ('W 


T - 

u.  ~ Xo 
R ^ Ro 


(S') 

"“I 

/ fl-u*- 

uT/  l~uT  . 


These  factors  will  not  change  the  form  of  any  of  the 
^cen  partial  differential  equations  and  will  only  come  into 
consideration  when  the  analysis  of  the  data  is  required. 
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Appendix  f 


The  Comp-uter  Programs 

The  One-dlmenslonal  Computer  Program 

The  listing  for  the  one-dimensional  model  consists  of 
a main  or  calling  program  and  four  major  subroutines.  The 
calling  program  serves  to  input  the  initial  conditions  to 
the  solving  routine,  to  input  step  sizes  and  accuracies  to 
the  solving  routine,  and  to  output  the  new  values  at  the 
time  step  established. 

The  variables  used  in  the  calling  routine  are  as 
follows ; 

AL  - Array  containing  relative  error  tolerances 
AW  - Scratch  array 
DS  - Size  of  time  step 
DU  - Size  of  spatial  step 
EO  - Constant  Eo 
F - Array  of  time  derivatives 
H - Step  size  estimate 

IB  - Negative  nxamber  means  no  size  estimate  given; 
zero  minus  size  estimate  given 
N - Number  of  spatial  points 
NA  - Total  niunber  of  time  steps  to  be  taken 
ND  - Total  number  of  equations 
P - Array  to  output  Eint  at  each  spatial  point 
PO  - Momentum  scale  factor 
RT  - Tolerance  for  solving  routine 


RT  - Tolerance  for  solving  routine 
S - Time  at  calculation 
SI  = Time  into  solving  routine,  initial  time 
SO  - Time  out  of  routine,  time  you  desire  answer 
SZ  - Array  containing  guesses  for  size  of  solutions 
U - Location  of  spatial  points  (scales) 
y - Array  containing  current  solutions 
YS  - Initial  input  and  final  output  array  listing 


P^D;?AH  QllTC<(I‘JPUT,0')TDUT) 

OIM-NSION  y S (17,4)  ,Y  (17,U,  12)  , - (17,4,12)  ,S(12)  ,AL  (68)  ,AW(68) 

l,S7(f3> , o (17) 

COMION  ^0,')U 

EXTZ2VJAL 

D'J  = . 1111111 

P.'‘A^q9,NA  , EC, OS, ^T, 30 

S0=-2v.  *N=17'5Nn=S9!n=-l 

CAUL  7E^o  (Y*;,l,‘10) 

ppIVT*,”ST4PTIM';  TIME=",S0,**INrERVAL=*V3S 

H=.:i*3S 

OOlt=l,NT 

1 AL(T)  = RT 
C02T=1, V 
U=(T-0. ) •*2/81, 

Y55(I  ,1)=-32.*(U**-^)4-12,*(U*U)-1,125*U+.015625 

2 CONTIMJE 

P?IMT98,iq,‘50,((Y3(I,J),J=l,4),I  = l,N) 
n07[ A=1,NA 
ST  = >0?S0=SD+0S 

cal.  q'.:<nO(PCM,'^'*,MO,SI,SO,Y';,5,Y,F,  AL,AW,H,SZ,I8) 
PPI'IT93,IB,SC(,  ((YS(I,  J)  ,J=1,4)  , 1=1,  N) 
no  3 1=1,17 
5 P(T)  =EEXT  (1,00) 

PPI'ITOS,  (P(I)  ,1=1,17) 

^ CONfTMUE 

99  F0P-(aT(r5,8F12.4) 

98  F0P8AT(7/l5,loE13,4/( 2X,1P4F13, 4)) 

95  F0P'^AT(2X,1PE13,4) 

ST03 

END 
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The  differential  equation  solving  routine  was  an 


eight-point  Newton-Rapsen  method  routine  and  will  not  be 
discussed  in  any  detail,  because  any  standard  differential 
equation  solving  routine  could  be  adapted  to  the  program. 

The  second  major  subroutine  consists  of  the  moment 
equations  derived  in  previous  sections.  This  subroutine, 
when  called  on  by  the  differential  equation  solving  routine, 
and  given  values  of  the  function  at  an  initial  time,  will 
return  values  of  the  time  derivative. 

The  following  is  a list  of  additional  variables  used 
in  this  subroutine: 

CS  - Hyperbolic  cosine  of  the  time 

DYi  - i = 1,2,3  derivatives  of  the  function  Yi , Ya, 

Y3  respectively 


RU  - A portion  of  the  scaling  parameter 

_ '^2 

Yi  - 1=  1,2,3  and  h values  of  N,  p,  P Eint  respect- 
ively. 


SU^’OJTINH  FCN(<?,y,F) 

DIMENSION  Y(17,^)  ,F(17,4) 

COMMON  FO,OU 
CS=-CDSH(S) 

CftL.  7£?0(P,17,4)  ^CftLL  7E^D  ("  (17,1)  , 17,  4) 

001T=1,17 

U=(T-D.)»*2/B1, 

RI)=50^T(  (l.-J)**3) 

R=C5*RJ 

Yl=y (t,l) ?y?=Y(I, 2)$Y3=Y(I,3)  *Y4=Y(I,4) 

0Y1=03DJ(Y,I,17)$DY2=090U(Y (1,2) ,1,17) S0Y3sDDDJ<Y (1,3), I, 17) 
F(T,  1)  (r?*DYl+DY2) 

F(I, ?)aO»  (Y3*DYl+DY3  + Y2*0Y2)fC3*(Y4  + EEXT(I,S) ) /EO 
FCI,  3)sR»  (2.*Y3*DY2<-Y2»0Y3) 

F(I,  4)  a-2. •F3*C'5* EXP  (Y1)*Y 2 

RETJPN 

ENO 
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The  third  subroutine  calculates  the  derivative  of  a 
function  by  an  eight-order  nmerical  routine.  The  n\amer- 
ical  differention,  as  described  in  Milne,  (Ref  20),  calcu- 
lates the  derivatives  at  eight  equally  spaced  points.  The 
accuracy  of  the  numerical  solutions  is  much  better  for 
points  in  the  center;  thus,  in  the  subroutine  the  coeffi- 
cients for  the  fifth  point  are  used  for  positions  5 through 
13  in  our  17  position  scheme.  A typical  derivative  at  the 
first  point  is  given  by 

t - (-llZloy,  + 1335-8  y,  + |283<?y.  + U04.V1 

-74/iy  - 475-0  -f  y>.  * ' 3’7V  ) 


where  h is  the  value  of  the  constant  interval. 

The  coefficients  to  the  numerical  solutions  are  stored 
in  the  array  RP(9,10),  with  the  last  member  of  each  array 
row  going  to  the  factor  which  multiplies  the  constant  in- 
terval size.  The  actual  calculation  of  the  derivative  is 
done  by  matrix  multiplication  of  the  proper  rows  of  Y,  with 
the  proper  rows  of  RP,  and  then  divided  by  RP(I,10)  an- 
other constant  interval  distance. 

The  following  new  variables  appear  in  DBDU: 

RP  - the  coefficients  as  discussed  above 

DBDU  - the  value  of  the  derivative  at  the  point  of 
interest 
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PUN^TIOM  0!!10U(Y,I,N) 

OIMrNSION  r(l),RP(q,l'?) 

COH'ION  "OtDU 

^ftT^  ( (?P(I,J)  ,J=1, 11)  ,1  = 1,9)  / 

1 -19  210.  , 133<5  6.  ,12839.  ,120  6., -7416.  ,-6750.  ,1633.  , 8314, , -■jg74.  , 1544 

2 4., -719  8., POC., 3 323. ,2844. ,124 2., -198., -785.  ,-422., 3 94.  ,15444.  , -8 

3 14. ,-4  4 0 8. ,-1615. , 2268. ,4^69. , 2736. ,-65  3. ,-29C  C. ,1318. , 15444. , 145 

4 4.,— V534»f— 3163. .493. ,30C6., 3024. ,841.,  — 1389. ,313.,! 5444 .,86.  ,—14 

5 2., -193. ,-126., 0., 126., 193.,  14 2. ,-9  6.,  1188., -310. ,13  98. ,-841.  ,-30 

6 24. ,-33  3 6. ,-43  0. ,3163. ,45  34.  ,-1454. , 15444. ,-l 319 ., 29 CO ., 653. , -27 3 
Y 6., -40  68., -2269.  ,1615., 44^8.  , 3 14.  , 15444. , -394. ,4  22. ,785. ,193., -12 

8 42.  ,-2  9 44.  ,-3323 ., -8r,  o ., 7 1 99 .,  15444.  , 3974.  ,-9314.  ,-1633.  ,675  0.  , 74 

9 16, ,-12  0 6. ,-12  93 P. ,-13359.  ,19210.  , 15444.  / 

IF(T-3)  1,2,5 

5 IF(v-r-4)  3,2,2 

1 nonjsVPgoKROd,!)  ,9,Y(1)  ,1  ,9,XX)/(OU»RP(I,10) ) 

GO  TO  4 

2 n5nj=y/PR0  0(PO(5,l)  ,9,Y(I-4)  ,1,9,  XX)  / (0U*RP(5,  10  ) > 

GO  '0  4 

3 090J  = \/t»?0  0(PP(I-8,l)  ,9,Y(9-S)  , 1 , 9,  XX ) / ( 9J*PP( I- 9, 10)  ) 

4 RiTJPN 
EMO 


The  last  subroutine  calculates  the  value  of  the 
electric  field  due  to  the  beam  at  a specific  point  within 
the  material.  Eq  (B-15)  is  the  source  for  values  returned. 
The  following  variables  are  found  in  EEXT: 

B = the  relativistic  constant 
C = the  speed  of  light 
E = the  value  of  the  electric  field 
EPSO  = to\  the  permittivity  of  free  space 
PI  = 3.14159 

Q = the  charge  in  one  pulse  of  the  beam 
R - radius  of  interaction  - a constant  in  1-D 
To  - time  scale  parameter 
V = velocity  of  beam 


139 


■ distance  scale  parameter 


fun: 

pi=? 
Q=i. 
C = 3. 
V=2. 

' ps: 
T0  = ’ 

xo=. 
R=.n 
9=V' 
U=  (I 
0 = 0/ 
A=V» 
r=n* 
+/ ( (( 
£Exr 
RETJ 
END 


TIOM 

.l'4l6 

r-os 

PE3 
939E3 
= 3.35 
• ?5  E- 

01 

-9)  /9 
('♦»’! 
T0»5:i 
((?»9 
2*3/  ( 


EEXT (I,X) 


4E-12 

17 


•EP5D) 

MH(X)-XO*U/SORT(1-U*U) 

/(l-3*3»  )-fl/S0RT(?»R*(  1-3*3) +A»A) ) 
1-3*3)  )*A-SaRT(R*R*(  l-0*3)+A*ft))**2) 


e Two-dimensional  Computer  Program 


The  components  of  the  2-D  program  are  identical  to  the 
1-D  model  except  for  the  extension  in  thej:*adial  direction 
of  the  spatial  nodes  and  the  electric  and  magnetic  fields 
The  extensions  require  additional  storage  space  and  sub- 
routines for  taking  derivatives  and  calculating  the  radial 
component  of  the  external  electric  field. 

The  main  program  serves  the  same  functions  as  described 
for  the  1-D  program.  The  following  is  a listing  of  addi- 
tional variables  used  in  the  2-D  main  program, 

DW  - The  step  size  in  the  radial  direction 
EE  - The  charge  on  an  electron 


EO  - Eo,  permittivity  of  free  space 


I 


RO  - The  radial  distance  parameter 

TO  - Time  parameter 

YY  - Conmion  storage  array 

Zo  - Axial  distance  parameter 

Z5  - The  unsealed  time  of  solution 

EM  - The  electron  mass 


P^C;9A*i  THiSIS  (IMPJT,Ot)T»UT) 

OTM-NSION  Y«;{17,q,10)  ,S(1?)  ,4L(  1513)  , AW  (1530)  , SZ  ( 1530  ) 
r-(1530,12)  ,YY(15ir,l?) 

RPn, 10)  ,‘=:0,r)U,n4,^3,70,EM,PE,TO 

EY’iP-^AL  P2N 

0U=?W=,11111 

N=1Y 

M0=1533 

I5=-l 

EE  = 1.6E:-1P 
P,0=o  » 25E-5 
70=. 27S5 
T0=7.25E-17 
EM=^,11E-11 

REA3  gBjNAjFOiOStPTtPO 

so=: 

CAL,  7ER0(Y«;,1,N1) 

H=,' 1*DS 
03  I 1=1, N3 

1 AL(D=^T 

REA3  1C,  ( (9P(I,  J)  ,J=1,1[))  ,T  = 1,5) 

10  POR'IA^dOFS,!) 

03  ? 1=1,17 
DO  2 J=l,3 
U=(T-3) »*2/5l, 

2 Y«5  (I  , J,l)  =-32.  »(U*»3)  +12.»(U’MJ)  -1,  125*U  + . L 15625 
PPI'IT  93,  ((  (Y^d,  J,<)  ,K  = 1,10)  ,J  = 1,9)  ,1=1,17) 

00  1 T=1,NA 

51=30 

S0=50f05 

CAL'.  3LCK33  (FCN,  NO,  WO , SI , 59 , YS,  5 , YY  , FF,  AL  , AW,  H,  SZ,T8) 
.PPIWT»,I9,S0 
75=T0*SINH(50) 

PRINT*, “TIME  = ’STS 

PRINT  93,  (( (YS(T, J,<) ,K  = 1,1C)  ,N  = 1,9)  ,L  = 1,17) 

1 CONTINUE 
99  Fo=»«iAT(2X,lP10E13.4) 


99  F0RNAT(T5,6E12.4) 

STO’ 

ENO 


The  differential  equation  solving  routine,  and  the 
derivative  subroutine  in  the  axial  direction  remain  the  same 
as  in  the  1-D  program. 

The  third  subroutine  which  calculates  the  values  of 
the  time  derivitives  serves  the  same  function  as  in  the  1-D 
program. 

The  following  variables  have  a different  meaning  or  are 
absent  in  the  1-D  model: 

B2-B3  - Hyperbolic  cosine  of  and  Pi^,  respectively 
CON  - The  constant  C2  described  in  Appendix  F 
EN  - Initial  number  of  electrons  in  the  plasma 
R - Value  of  the  radius  (scales) 

RW  - A portion  of  the  derivative  with  respect  to  R 
W - Location  of  nodes  in  the  scales  R direction 
X - Time  (scaled) 


SUn^OJTIME  FCN  (X,Y,c-) 

0I'<-N5I0N  Y (17 ,9,  1C»  , F(17, 9»10)  , FORCE  (17, 9 » 2)  ,090W(8,  IG  ) 
♦ ,0P'nj(17, 1C) 

COMMON/OME/  RPO,lC)  , -0,0U,  O'!  ,R0,Z9,  FM,EE,TO 
NW  = 5 
N'J  = 17 

CON=2.4BOE7 
C=7.  E9 
EN=1,9E4-23 
no  1 J=2,9 
00  1 1=1,17 
U=(r-9,>  /9. 

W= (J-1,) /9. 

«?=OT*M/(  1,-H) 

RU=30RT(  (1-U*U)»*?)  fir) 

RW=(1-W)**2/R0 
ft  = n«?H(X) 

Y1=Y  (I,  J,  1) 

Y2=Y  (I,  J,2) 

Y3=Y  (T,J,3) 

YUrY(T,J,4) 

Y5=Y (I, J,5) 

Y5=Y (I, J,6) 

Y7=Y  (I,  J,7) 
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Y8='f  (I,  J,8) 

V9  = Y <r,  J,9) 

Y1  1= Y(I, J,10) 
f^?  = 'n9H(Y  (T,  J,2) ) 

93=:0?H(  Y (T  , J,3)  ) 

(J,l)='*c;9H(Y(T,l,  1)  , 

(j,  ?)  =r>rPw(Y(  T,  1,  ?)  , J,MW) 

D9n<*  (J,3)  =T'-9W(Y(  1, 1, 3)  , J,  NW»  •oW/'»3 
O'^OJ  (T,l)=T3;oU(v(l,J,l)  ,I,NU)*9'J*R3r‘^0 
ORf'J  (1,3)  =9rRU(Y(  1,  J,  ?)  ,1,  *Rl)*R9/  (Z9*3?) 

nmi  (t ♦ 3)  =orou(Y(  1,  j,  •»)  ,i,  *ju>  *R'j*ro/  C93»zo) 
no  ’ Y=u,ic 

(J,O=0rRW(Y(T,l,K)  , J,NW)  »=»M 
? n9nj(T,<)=nrR'j(Y(i,j,'<)  *rj*ro/zo 

FOO^Fd,  J,l)  = -(  Y8-V»*Yi:.)  •COM 
Fo^crd  , J,Z)  =-C0M*  (V9  + V2*Y1C) 

T=Y?*09n'M  ( J,  1>+-Y3*n=«n'l(l,i)  (J,2)  *0900  (1,3)  ♦Y2/R 

Fd,  J,l)  =a*  (Oaoi-H  J,l)  •Y3/Y1+1‘^0  Jd,l)  •Y3/Y1  + Y1*(0B0M(  J,2)  ♦0900(1,3 

♦ )^Y?^R)  ) 

Fd,  J,2)  =ft»FnRrC  d,  J,  l)-6*  ( YV*030W(  J,  1)  +Y5»  0900  (1,1)  ♦090H(  J,4) 

♦ ♦0930(1 ,5)  ♦Yt/R-R»YP.  + Y2^090W(J,  2) +Y  3»090  J ( 1 , 2 > > 

F(I,  J,  3)  =«»FO9r'd,J,2)-a*(Y9“O90W(J,l)^Y7»09O'J(I,l)^090M(J»5) 

♦ ♦0®3U(I,7)+Y?*n90'.'l  (J,  (1,3)  ) 

P(I,  J,  4)  =-4»(T^Y9*Y2+2.*(n9nM ( J , 2) ♦ 0900 ( 1 , 3 ) ) * ( Y4^Y5^ 

♦ Y?»Y?fY4»Y2)+Y2*0  90M  ( ),£*>  ♦ Y 3*030U ( I , 5) ) 

F(I, J,5) =-a*(T»Y2*Y3+ (Y5  + Y7  + Y2» Y3^Y3*Y3)  • (090W(  J,  2)^ 

♦ 0«OJ  (!,?))♦  (Y4^Y9+Y2»  Y2  + Y2*Y3)  • (093W  ( J,  0)  4-093lJ(I,  3) ) ♦ Y2»090M  ( J, 5 > 
♦♦Y3»09DJ (1,9) ) 

F(  I,  J,5)  =-&♦(  YP.^YZ/R  + YZ^OBf'M  ( J,fi)  ♦4,*Y6/2)  ♦Yl* 030'J ( I ,6)  ) 

Fd,  J,7)  =-ft*  (T*Y3*Y3^2»  ( Y5+ V7  + Y ?•  Y3  + Y 3»  Y3)  • (090W  ( J,  3) 

♦ ♦ORTUd  , 3))  + Y2*09  0W(  J,7)+Y3»090J  (T,  7)  ) 

Fd,  J,  ft)  =A»  (-0900(1,1.  ) ♦Y2^EYP(  Yl)  ) 

Fd,  J,9)  = a*  (099M(  J,10)^Y1?'=>>Y3*EXP(Y1)  ) 

FCN  74/74  0PT=1  FIN  4, 5^414  09/ 

F(I,  J,1C) =A*(0O0W(J,9> -0900  (1,3) ) 

1 CONTTNOE 
RETIRM 
END 


The  subroutine  DERW,  takes  the  derivative  of  the  func- 
tion in  the  radial  direction.  This  numerical  derivative 
is  taken  in  exactly  the  same  manner  as  was  the  radial  com- 
ponent, except  the  matrix  multiplication  is  done  along 
the  radial  direction  of  the  Y matrix. 
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F'JM:tIO'‘J  nERX(Y,J,N) 

OT'ITN^rn*^  VC^) 

COM'^OM/ONE/  RPCT.IO)  , EO,  lU,  , ZO  , ZM , FE , TO 

0£R-J  = \/PR0n(RP(  J,l)  ,9,Y(1)  , 14'*,9,XX)  / (OW*RP(  J,13)  ) 
pryjosi 
E^n 


The  external  electric  field  in  the  axial  direction  is 
calculated  in  a similar  manner  as  was  done  in  the  1-D  pro- 
gram except  the  radius  of  the  point  is  allowed  to  vary. 

The  external  electric  field  in  the  radial  direction  is 
obtained  through  the  use  of  Eq(S-l4).  The  variables  in 
EEXTR  are  consistent  with  the  rest  of  the  program. 


PUrnriON  ZEXTR(T,J,X) 

OOf1'<OM/ONE/  RP{9,10)  , "0 , 00,  , R D ,ZD , EM, EE , TO 

OCTS  PI,a,:,Y,EPS0/3, 1U159, 1. 0- -6,3. 0Z8 , 2. 999E8 , 8 . 854E-12/ 

9 = V''' 

W=(J-l,)/9, 

P.=  PD*W/(1,-W) 

U=(T-9,)/9. 

Z=Z9*J/SQP.T(t.-M*|j) 

0=0/  (4*PI*s’P<;0) 

H = 7»TO*^INM(X)-TO*U/SOPT(1-U*  IJ) 

F=S0PT(9»P*  (l,-9*9)  ♦A^A) 

G=?. *9/ (1,-P*P) 

E=0*P**(l,-»»*9)/(  ((G*A-P)**2)  *F) 

EFxrPs? 

P‘^TJPM 

EMO 


^kk 


Appendix  G 


Derivation  of  Plasma  Frequency  From  The 
Moments  of  the  Vlasov  Equation 

Because  we  have  chosen  to  use  the  plasma  oscillation 
as  a standard  for  validation  of  the  technique  and  models, 
it  becomes  necessary  to  justify  that  u)p  (the  plasma  fre- 
quency) is  a valid  solution. 

If  we  linearize  the  first  three  moments  of  the  Vlasov 
equation,  Eqs  (D-6),  (d-13)  and  D-l?)» 

N = No  + Ni 
P = Po  + Pi 
F = Fo  = Fi 

f\j  '\,2  ^2 

P2  = Po  + Pi 

-^3 

P = 0 

where  No,  Po,  Po  Fo  are  constants  and  Ni,  P’1,  Pi,  Fi  are 

small  perturbations  to  the  initial  values,  and  if  we  as- 

'^2 

sume  the  plasma  to  initially  be  cold  (Po  = o,  P = o and 
Fi  = o;  no  Initial  momentum  or  temperature)  the  first  three 
moments  become 
First  Moment 


Second  Moment 


4-  ( A/o  + A/,)  R ^ p,  ^ 

Third  Moment 

r*  +//,) 0"  ^ I ^ 

(c  0 

If  we  now  keep  only  the  first-order  perturbation  terms 
and  assume  that  a small  perturbation  times  a small  pertur- 
bation is  negligible,  Eqs  (1),  (2)  and  (3)  become 

- --k  kr 

-V*  ^ ~ A(^‘  ir 


A/o  m O. 

From  Eq  (6),  we  know  that  No  r(:o  thus,  Xf  ^ 
or  perturbations  oftemperature  in  a plasma  are  not  a function 
of  time.  Eq  (5)  then  reduces  to 

M,  ^ yV.  . (G  - l) 

Let  us  now  look  at  the  force  term.  If  we  assume 
Ni,  Pi,  Fi  have  the  form 
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-i(kx-ut) 


then 


where 


But  we  know 


Ai  = Ai  e 
F = e(E^  + E‘*) 

E^  = the  external  electric  field 
E^  = the  internal  electric  field. 


T 


e 


^ I’  --I^u 


Pi  a/ e 


and  that 


All 


~ C Cj 


Substituting  this  result  into  Eq  (7)  yields 


No 


C A/a  ^ ^ £*//• 

A«c  < 


{'g  '?) 


Taking  the  partial  derivative  of  Eq  (q)  and  canceling  like 
terms  yields 


p -_e_E 


uJai<  f • / 

Substituting  this  equation  into  Eq  (4)  and  reducing,  de- 
fines Ni  as 


‘ - A4»C*) 


Thus,  if  there  are  no  density  perturbations  involved  in 

simple  harmonic  motion  N,  = o and  Ex  = o 

u)  = Noe*  = oj*p 
Meo 

where  u/p  = the  plasma  frequency.  Ref.  15  demonstrated  the 
fact  that  as  the  plasma  becomes  finite  radially,  instead 
of  the  infinite  model  chosen  in  the  deviation,  the  frequency 
of  the  plasma  approaches 

UJ  - ^ 

rr  ' 
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